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ABSTRACT 

Gravitational-wave (GW) recoil of merging supermassive black holes (SMBHs) may influ- 
ence the co-evolution of SMBHs and their host galaxies. We examine this possibility using 
SPH/N-body simulations of gaseous galaxy mergers, including SMBH accretion, in which 
the merged BH receives a recoil kick. This enables us to follow the trajectories and accretion 
of recoiling BHs in self-consistent, evolving merger remnants. In contrast to recent studies 
on similar topics, we conduct a large parameter study, generating a suite of over 200 sim- 
ulations with more than 60 merger models and a range of recoil velocities (v^). With this, 
we can identify systematic trends in the behavior of recoiling BHs. Our main results are as 
follows. (1) BHs kicked at nearly the central escape speed (v csc ) may oscillate on large or- 
bits for up to a Hubble time, but in gas-rich mergers, BHs kicked with up to ~ 0.7 v csc may 
be confined to the central few kpc of the galaxy, owing to gas drag and steep central poten- 
tials. (2) u csc in gas-rich mergers may increase rapidly during final coalescence, in which 
case recoil trajectories may depend sensitively on the timing of the BH merger relative to the 
formation of the central potential well. Delays of even a few x 10 7 yr in the BH merger time 
may substantially suppress recoiling BH motion for a fixed kick speed. (3) Recoil events 
generally reduce the lifetimes of bright active galactic nuclei (AGN), but short -period recoil 
oscillations may actually extend AGN lifetimes at lower luminosities. (4) Bondi-Hoyle ac- 
cretion dominates recoiling BH growth for v-klv csc < 0.6 — 0.8; for higher v^, the BH accretes 
primarily from its ejected gas disk, which we model as a time-dependent accretion disk. (5) 
Kinematic ally-offset AGN may be observable either immediately after a high-velocity recoil 
or during pericentric passages through a gas-rich remnant. In either case, AGN lifetimes may 
be up to ~ 10 Myr (for Av > 800 km s _1 ), though the latter generally have higher luminosi- 
ties. (6) Spatially-offset AGN can occur for Wk ;^0.5 — 0-7 v csc (for AR > 1 kpc); these 
generally have low luminosities and lifetimes of ~ 1 — 100 Myr. (7) Rapidly-recoiling BHs 
may be up to ~ 5 times less massive than their stationary counterparts. These mass deficits 
lower the normalization of the Mbh — a * relation and contribute to both intrinsic and overall 
scatter. (8) Finally, the displacement of AGN feedback after a recoil event enhances central 
star formation rates in the merger remnant, thereby extending the starburst phase of the merger 
and creating a denser, more massive stellar cusp. 

Key words: black hole physics - gravitational waves - accretion, accretion disks - galaxies: 
kinematics and dynamics - galaxies: evolution - galaxies: active 



1 INTRODUCTION 

Within the hierarchical structure formation paradigm, a significant 
fraction of galaxy growth occurs via successive mergers. Given the 
ample evidence that most, if not all, local galaxies host central 
supermassive black holes ( SMBHs; e.g., [Kormendy & Richstone| 



[T995l[Richstone et al.| 19981 |Ferrarese & Ford|2005| >, major galaxy 
mergers (those with mass ratios q > 0.25) should inevitably result 
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in the formation of SMBH binaries. The fate of these binaries is 
somewhat uncertain and likely depends on the conditions within 
their host galaxies. In highly spheroidal, gas-poor galaxies, these 
binaries may "stall" at separations of about a parsec for up to a 
Hubble time (e.g., |Begelman et al.|1980||Miiosavljevic & Merritt| 
2001 ; Yu 20021. In gas-rich galaxy mergers, however, gas is driven 
to the central region of the merging galaxies simultaneously with 
the formation of the SMBH binary. Numerical simulations indi- 
cate that BH merger timescales may be much shorter than galaxy 
merger timescales in this case (~ 10 6 — 10 7 yr from the hard binary 
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stage; e.g.,;Escala et al 2005j |Dotti et al.|2007| l. This implies that 
the SMBH binaries most able to accrete gas and produce electro- 
magnetic signatures are likely to be short-lived. 

To date, observations seem to confirm this scenario. Several 
quasar pairs with kpc-scale separations have been found in merg- 
ing galaxies ( Komossa et al. 2003 ; Bianchi et al. 2008 , Comerford 
|et al.|2009b| |Green et al.|2010f 7Recently, spectroscopic and pho- 
tometric observations have found that between 10 ~ 3 — 10 -2 of 
active galactic nuclei (AGN) at z < 0.7 may in fact contain dual 
BHs with separations of ~ k pc {Comerford et" al. 2009a; ^Smith| 
|et al.|2010|[Liu et al.|2010blal|Smith et al.|2010||Fu et al.|2010| >. 
However, only one unambiguous case of a SMBH binary has been 
found, with a separation of ~ 7 pc (Rodriguez et al. 2006). Several 
additional objects have been identified as candidate SMBH bina- 
ries but are still unconfirmed (Sillanpaa et al. 1988, Komossa et al 



2008 



|Dotti et al. ;2009; Bogdano vTc et al.|2009f |Boroson & Lauer 



2009). Thus, it is clear that the large majority of AGN do not con 
tain binary SMBHs, which suggests that any substantial population 
of long-lived SMBH binaries must exist in gas-poor environments 
where they are quiescent. 

If most binaries are therefore assumed to merge on reason- 
ably short timescales, then gravitational-wave (GW) recoil must be 
a common phenomenon throughout the merger history of galaxies. 
GW recoil is a natural consequence of gravitational radiation from 
BH binary mergers jPeres|1962)|Bekenstein|1973|ptchett & De-| 
|tweiler|1984| >. If the binary system has any asymmetry - unequal 
masses, spins or spin orientations - then gravitational waves are 
radiated asymmetrically, resulting in a net linear momentum flux 
from the final BH at the time of merger. This causes the BH to 
recoil in the opposite direction. Whether these recoil kicks were 
large enough to be astrophysically interesting was uncertain until a 
few years ago, because accurate calculations of the recoil velocity 
require BH merger simulations using full general relativity. These 
simulations have revealed that GW recoils may be quite large. Re- 
coil velocities up to 4000 km s~ are possible for special configura- 
tions, which is far greater than galactic escape speeds (Campanelli 
|et al.|2007a|b| >. The implication that SMBHs may spend substantial 
time in motion and off-center has opened a new line of inquiry into 
the ramifications for SMBHs and their host galaxies. 

A useful starting point for such inquiries would be the distribu- 
tion of recoil kick velocities, but this is quite difficult to ascertain in 
practice. The final velocity depends sensitively on not just the mass 
ratio of the progenitor BHs, but on their spin magnitudes and orien- 
tations as well. The distributions of SMBH binary mass ratios and 
spins at various redshifts have been estimated using halo merger 
trees and semi-analytical models o f SMBH growth (e.g., | yolonteri| 
|et al.|2003|[2005||King et al.|2008|[Berti & Volonteri|2008> . These 
distributions depend on a number of model assumptions, however, 
and the BH spin orientations prior to merger are far more uncertain 
still. Thus, based on recent results from BH merger simulations us- 
ing full numerical relativity (NR), several groups have calculated 
kick probability distributions as a function of BH mass ratio for ei- 
ther fixed or random values of BH spin, with the assumption that 
the spins are randomly oriented jSchnittman & Buonanno||2007 



Campanelli et al.|2007a||Baker et al.|2008[|Lousto et al.|2010a]b 



van Meter et al.|2010| >. Their results are in good agreement with 
each other and imply that the fraction of high-velocity GW recoils 
is substantial. This underscores the potential for recoil events to be 
an important component of galaxy mergers. 

It is quite possible, of course, that the spins of SMBH binaries 
are not randomly oriented, but are preferentially aligned in some 
way.jBogdanovic et al. (2007) have suggested that torques in a cir- 



cumbinary gas disk may align the BH spins with the orbital axis 
of the disk. In this case, the resulting in-plane kicks would have 
a maximum recoil velocity of < 200 km s _1 , although spins that 
became anti-aligned by the same mechanism could result in recoil 
velocities up to 500 km s _1 (e.g., Gonzalez et al. 2007 Campanelli 



et al. 



ct al. 



2007a; |Baker et al.|2008) . Additionally, simulations by |Dotti 



(2010) demonstrate efficient spin alignment of merging BHs 
in the presence of a highly coherent accretion flow, although it is 
unclear how efficient this process might be in a gaseous environ- 
ment that includes, e.g., star formation. jKesden et al.|p010) have 
recently demonstrated that BH spin alignment may instead occur 
via relativistic spin precession, regardless of whether a gas disk is 
present. The aforementioned recoil kick probability distributions 
are therefore upper limits on the actual distributions. 

Numerous possible consequences of GW recoil events have 
been discussed in the literature. GW recoils may have a large effect 
at high redshifts, where escape velocities of galaxies are smaller 
(e.g., |Memtt et al.|2004[|Madau & Quataert|2004[|Volonteri|2007l >. 
This is a concern for attempts to understand, for example, the origin 
of the 2 = 6 SDSS quasars (e.g., |Fan et al.|2001|[2053) . |Volonteri| 
|& Ree s ( 2006 ) suggest that growth of SMBHs at high z must oc- 
cur only in highly-biased halos. Using cosmological hydrodynamic 
simulations, Sijacki et al. ( 2009 1 investigate BH growth in massive 
high-z halos, including the ejection of BHs with recoil velocities 
above v csc . Their findings are consistent with the observed popula- 
tions of bright quasars at z = 6, despite the effects of GW recoil. 

At lower redshifts, recoiling BHs may produce electromag- 
netic signatures. The main signatures we will focus on here are 
recoiling AGN that are either spatially or kinematically offset from 
their host galaxies (Madau & Quataert 2004, Loeb 200 7j|Blecha &| 
|Loeb|2008||Komossa & Merritt|2008a||Guedes et al.|2009) . Other 
possible signatures include flares from shocks induced by fallback 
of gas marginally bound to the ejected BH (|Lippai et al.||2008[ 
|Shields & Bonning|[2008l |Schnittman & Krolik||2008|>, enh anced 
rates of stellar tidal disruptions ( Komossa & Merritt 2008b j [Stone| 
|& Loeb||2010 i, and compact stellar clusters around ejected BHs 
( jO'Leary & Loeb|2009||Merritt et al.|2009) . 

Thus far, no confirmed GW recoil events have been observed. 
An inherent challenge in observing offset quasars is that larger spa- 
tial or kinematic offsets are easier to resolve, but less gas will be 
bound to the recoiling BH at higher recoil velocities, so its AGN 
lifetime will be shorter. Bonning et al. (2007) conducted a search 
for kinematic offsets in SDSS quasar spectra and found a null result 
at a limit of 800 km s _1 . Several recoil candidates have been pro- 
posed, but their extreme inferred velocities should be exceedingly 
rare. Indeed, the recoil candidate with a 2600 km s _1 offset, pro- 
posed by Komossa et al. (2008 ), may in fact be a superposition of 
two galaxies ( Heckma n et al.|2009|| Shields et al.|2009a) o r a binary 
SMBH system ( |Dotti et al.|2009||Bogdanovic et al.|2009| >. Another 
candidate with an even higher (3500 km s~ ) offset is most likely 
a double-peaked emitter (Shields et al. 2009b). Recently, |Civano| 
|et al.| f2010) have suggested that an unusual galaxy discovered in 
the COSMOS survey by |Comerford et al.|p"009b) , which those au- 
thors proposed to be a dual SMBH system, may in fact be a re- 
coiling BH, as new spectra indicate a kinematic offset of 1200 km 
s _1 . This candidate has a less extreme velocity than the others, but 
further observations are needed to confirm the nature of this object. 
Additionally, the SMBH in M87 has recently been observed to be 
spatially offset by ~ 7 pc, which could possibly be explained by a 
past recoil event (Batcheldor et al. 2010}. 

In addition to producing direct observational signatures, GW 
recoil may play a role in the co-evolution of SMBHs and their 
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host galaxies. Strong empirical correlations exist between SMBH 
mass and properties of the host galaxy bulge, including the bulge 
luminosity, mass, and stellar velocity dispersion (e.g., Kormendy 
1995] pagbrrian et al.|| 19981 |Gebhardt et al.||2000| 



& Richstone 



Ferrarese & Merritt 2000; Merritt & Ferrarese 2001, Tremaine 



et al.| [2002; Marconi & Hunt 2003). These correlations are well- 



reproduced by galaxy formation models in which much of BH and 
galaxy bulge growth occurs via successive mergers, and in which 
merger-triggered BH fueling is self-regulated via AGN feedback 



(e.g., |Silk & Rees][T998"l |Wyithe & Lo eb 2003 ; Pi Matteo et"aT] 



[20051 |Hopkins et al.|2006a| T However, because GW recoil events 
may occur simultaneously with this rapid BH accretion phase, re- 
coils could significantly disrupt the coordinated growth of BHs and 
galaxy bulges. In particular, GW recoil may contribute to scatter in 
the Mbh — cr* relation caused by ejected (|Volonteri|20071 > or bound 
recoiling < |Blecha & Loeb|2008||Sijacki et al.|2010| > BHs; we will 
examine the latter possibility in greater detail. It is also unclear a 
priori what effects, if any, GW recoil may have on the host galaxies 
themselves. In purely collisionless galaxies, Boylan-Kolchin et al. 
( 2004 1 and |Gualandris & Me rritt ( 2008} have shown that bound re- 
coiling BHs may scour out a stellar core. In gas-rich galaxy merg- 
ers, copious of amounts of cold gas are driven to the central galactic 
region during coalescence, triggering a luminous starburst such that 
the galaxy may appear as a ULIRG (Sanders et al. 1988a; Sanders 
& Mirabel 1996). Feedback from a central AGN may terminate 
the starburst phase by expelling gas and dust from the central re- 
gion (e.g., |Hopkins et aL|2006a| |2008b[ |Somerville e"tluT|2008] >. 
We will investigate whether observable starburst properties or gas 
and stellar dynamics may be affected by the sudden displacement 
of this central AGN via recoils (see also recent work by Sijacki 
|et al.|2010| >. 

Our current work was initiated as a follow-up of Bl echa &] 
|Loeb| {2008) >, hereafter BL08. BL08 explored the trajectories and 
accretion of BHs on bound orbits in a static potential with stellar 
bulge and gaseous disk components. They found that recoil kicks 
< n csc could produce long-lived ( > Gyr) oscillations of the BH; 



similar results were found by Gualandris & Merritt ( 2008 1, Kornre- 
|ich & Lovelace] p008l > & |Guedes et al.| ( |2009F . BL08 also demon- 
strated that GW recoil can affect SMBH growth even when the BH 
is not ejected entirely from the host galaxy. By removing accreting 
BHs from the central dense region and thereby limiting the BH's 
fuel supply, BHs receiving recoil kicks > 0.5 v esc may be less mas- 
sive than stationary BHs. 

In the present study, we use hydrodynamic simulations to self- 
consistently calculate the dynamics and accretion of recoiling BHs 
in realistic merger remnant potentials. We simulate galaxy merg- 
ers using GADGET- 3 ( |Spring el 2005 1, an SPH/N-body code, and 
apply a recoil kick to the BH at the time of BH merger. Due to 
the detailed initial conditions and many free parameters involved in 
such merger simulations, we conduct a large parameter study with 
dozens of galaxy merger models and a wide range of kick veloci- 
ties. Each merger model is simulated with at least one kick velocity 
and also with no recoil kick, for comparison. With this suite of re- 
coil simulations, we are able to observe trends in the behavior of re- 
coiling BHs in different environments. When paired with our time- 
dependent, sub-resolution models for recoiling BH accretion, this 
approach allows us to, for example, determine velocity-dependent 
recoiling AGN lifetimes and estimate the effect of GW recoil on 
scatter in the A/bh — relation. 

As this paper was in the final stages of preparation, two papers 
appeared that also involve hydrodynamic simulations of recoiling 
BHs in galaxy merger remnants (Guedes et al. 2010, Sijacki ~t al.| 



|2010[ >. However, each study used only a few simulations, probing 
only a small part of parameter space. Guedes et al. (2010) used re- 
sults of three merger simulations as initial conditions, and ran hy- 
drodynamic simulations of recoiling BHs in these merger remnants 
for a short time. These initial trajectories were used to calibrate 
semi-analytic calculations of the recoil trajectory in the remnant 
potential. A range of kick velocities was tested in each of the three 
merger models used. Sijacki et al. (2010) simulate recoils in an iso- 
lated, stable disk galaxy; they also use this galaxy model as the 
initial condition for a full merger simulation with GW recoils. Both 
studies also include prescriptions for accretion onto recoiling BHs. 
We provide a comparison of our results to the findings of these pa- 
pers in § [6] 

Our simulation methods, GW recoil treatment, and galaxy 
merger models are outlined in § [2] § [5] - [5] contain our results, 
organized as follows. In § [3j the general characteristics of our 
galaxy merger simulations with stationary central BHs (i.e., no re- 
coil kicks). We discuss the variation in merger dynamics and rem- 
nant morphologies between models. In §j4] the dynamics of recoil- 
ing BHs are discussed. We describe the general characteristics of 
recoil trajectories and trends between models in § |4.1| In § |4.2| we 
analyze the sensitivity of recoil trajectories to the BH merger time, 
and in § |4.3| we investigate the dependence of recoil trajectories 
on the direction of the kick. §[3] contains our results for the accre- 
tion and feedback of recoiling BHs. In § |5.1| we examine accretion 
and feedback of recoiling BHs, compare to the fueling of stationary 
BHs, and calculate AGN lifetimes from Bondi-Hoyle accretion. In 
§ |5.2| we describe our time-dependent analytic model for calculat- 
ing accretion onto a recoiling BH from a disk of gas ejected with the 
BH. Using this accretion model, along with accretion from ambi- 
ent gas, we present in § |5.3| the velocity-dependent AGN luminosi- 
ties and active lifetimes for recoiling BHs. In particular, we calcu- 
late lifetimes for spatially-offset and kinematically-offset recoiling 
AGN. In § |5.4| we compare the Mbh — c* relations derived from 
our set of no-recoil simulations and our set of high- velocity recoil 
simulations, and discuss implications for the observed A/bh — cr* 
relation. Finally, in § |5.5| we explore the effects of GW recoil on 
star formation rates and the host galaxy structure. We summarize 
and discuss our results in §|6] 



2 SIMULATIONS 

For our simulations of galaxy mergers, we employ the smoothed 
particle dynamics (SPH) code GAD GET- 3 l|Springel|2005|>, which 
conserves both energy and entropy (Springel & Hernquist 2002). 
The code includes radiative cooling as well as a subresolution 
model for a multiphase interstellar medium (ISM) (Spri ngel &| 
Hernquist 2003) that accounts for star formation and supernova 
feedback. In addition, the code models BHs as gravitational "sink" 
particles that contain a BH seed and a gas reservoir. The reservoir 
is replenished by stochastic accretion of neighboring gas particles, 
but the actual accretion rate onto the BH is calculated smoothly us- 
ing the Bondi-Hoyle-Lyttleton formula with locally-averaged val- 
ues for the density and sound speed. Because the gas around the BH 
cannot be resolved at higher densities below the spatial resolution 
limit, the accretion rate calculated from these values is multiplied 
by a constant factor. Following other authors, we assume a value of 
100 for this factor (e.g. |Springel et al.|2005b||Hopkins et al.|2006b| 
|Johansson et al.||2009| l. Angular momentum is conserved during 
accretion of gas particles, but because this is a stochastic process 
owing to the finite mass resolution, we also introduce an analytic 
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Initial model parameters Resulting merger quantities 
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1024 


1112 


8.6 


q0.5fg0.4d 


0.5 


204.1 


8.368 


0.4 


7.1 


150 





-30 


45 


1.90 


1061 


1117 


5.2 


q0.5fg0.3a 


0.5 


204.1 


8.368 


0.3 


7.1 


30 


60 


-30 


45 


1.97 


993 


994 


0.1 


qO.5fgO.3b 


0.5 


204.1 


8.368 


0.3 


14.3 


30 


60 


-30 


45 


2.45 


920 


1003 


9.0 


qO.5fgO.3c 


0.5 


204.1 


8.368 


0.3 


7.1 


150 


60 


-30 


45 


1.97 


1012 


1017 


0.5 


q0.5fg0.3d 


0.5 


204.1 


8.368 


0.3 


7.1 


150 





-30 


45 


1.92 


1014 


1069 


5.4 


q0.5fg0.3e 


0.5 


204.1 


8.368 


0.3 


7.1 


90 


60 


-30 


45 


2.05 


915 


1126 


23.1 


qO.5fgO.3f 


0.5 


204.1 


8.368 


0.3 


7.1 


90 











1.94 


978 


1029 


5.2 


qO.5fgO.3g 


0.5 


204.1 


8.368 


0.3 


7.1 


60 


60 


150 





2.23 


805 


820 


1.9 


qO.5fgO.3h 


0.5 


204.1 


8.368 


0.3 


7.1 














1.84 


927 


1111 


19.8 


q().5fg0.3i 


0.5 


204.1 


8.368 


0.3 


7.1 


180 











1.84 


1028 


1146 


11.4 


qO.5fgO.3j 


0.5 


204.1 


8.368 


0.3 


7.1 


180 





180 





2.09 


823 


930 


13.0 


q().5fg0.3k 


0.5 


204.1 


8.368 


0.3 


7.1 


10 





-10 





1.84 


1046 


1087 


4.0 


q0.5fg0.3M10x 


0.5 


2041 


83.68 


0.3 


15.4 


30 


60 


-30 


45 


2.02 


2081 


2070 


-0.5 


q().5fg0.2a 


0.5 


204.1 


8.368 


0.2 


7.1 


30 


60 


-30 


45 


2.02 


830 


862 


3.8 


qO.5fgO.la 


0.5 


204.1 


8.368 


0.1 


7.1 


30 


60 


-30 


45 


2.12 


777 


798 


2.7 


q().5fg0.1b 


0.5 


204.1 


8.368 


0.1 


14.3 


30 


60 


-30 


45 


2.56 


812 


813 


0.1 


q0.5fg0.1c 


0.5 


204.1 


8.368 


0.1 


7.1 


150 


60 


-30 


45 


2.12 


811 


813 


0.2 


q().5fg0.1d 


0.5 


204.1 


8.368 


0.1 


7.1 


150 





-30 


45 


2.12 


814 


810 


-0.5 


q0.5fg0.1M10x 


0.5 


2041 


83.68 


0.1 


15.4 


30 


60 


-30 


45 


2.07 


1737 


1737 


0.0 


q0.5fg0.04a 


0.5 


204.1 


8.368 


0.04 


7.1 


30 


60 


-30 


45 


2.18 


740 


748 


1.0 


q0.5fg0a 


0.5 


204.1 


8.368 





7.1 


30 


60 


-30 


45 


2.30 


689 


706 


2.4 


q0.333fg0.4a 


0.333 


181.4 


7.438 


0.4 


7.1 


30 


60 


-30 


45 


2.28 


810 


816 


0.7 


qO.333fgO.3a 


0.333 


181.4 


7.438 


0.3 


7.1 


30 


60 


-30 


45 


2.38 


856 


858 


0.3 


qO.333fgO.la 


0.333 


181.4 


7.438 


0.1 


7.1 


30 


60 


-30 


45 


2.57 


741 


744 


0.3 


q0.25fg().4a 


0.25 


170.1 


6.974 


0.4 


7.1 


30 


60 


-30 


45 


2.94 


829 


830 


0.1 


q0.25fg0.3a 


0.25 


170.1 


6.974 


0.3 


7.1 


30 


60 


-30 


45 


2.90 


852 


852 


-0.1 


qO.25fgO.la 


0.25 


170.1 


6.974 


0.1 


7.1 


30 


60 


-30 


45 


3.18 


722 


725 


0.3 



Table 1. Galaxy merger models. Boldface entries denote those for which we have varied v k . Fiducial-mass merger models are labeled as q[value]fg[value][orb], 
where "q" is the galaxy mass ratio, "fg" is the initial gas fraction, and each letter orb corresponds to a specific orbit ("a" being our "fiducial" configuration). 
High- and low-mass models are denoted by q[value]fg[value]M\factor][orb], where "M" is the ratio of the primary galaxy mass to the fiducial mass. 
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accretion drag force (oc M v) calculated from the Bondi accretion 
rate at each timestep. These prescriptions are described in more de- 
tail in |Springel et"aT] i |2005bF . 

We note that, as in all hydrodynamic simulations involving BH 
accretion, the exact nature of the accretion on sub-resolution scales 
is unknown, so it is necessary to assume a sub-resolution model. 
While the validity of the Bondi-Hoyle-Lyttleton formula (with a 
multiplicative factor) in our framework is indeed an assumption, 
few constraints exist on the gas flow from large scales down to 
the BH. The Bondi-Hoyle-Lyttleton model provides a physically- 
motivated approximation of the gas captured by a BH, including a 
BH in motion, and its extensive use in the literature helps to place 
our work in a larger context of galaxy evolution studies. We further 
assume that the accretion factor remains constant even after a GW 
recoil event. This latter assumption is justified in that the spatial res- 
olution limit, which motivates the use of this factor, always remains 
constant. Furthermore, owing to the sink-particle treatment of BHs 
in the code, the accretion factor does not determine the amount of 
gas bound to the BH at each timestep; it influences only how this 
gas is accreted onto the BH. From a practical standpoint, using the 
same accretion prescription throughout enables a more direct com- 
parison between recoiling and stationary BHs - one of the main 
goals of this study. 

Because we are presently interested in GW recoil, we al- 
low for an arbitrary velocity to be added to the remnant BH at 
the time of the BH merger, i mrg . To more easily compare recoil 
events in different merger models, we scale this velocity to the 



fine Vesc = \] — 2$(xbh) at each timestep, i.e., the escape speed at 
the position of the BH, which is close to the position of the potential 

.|Campanelli et al.H2007a|>, 



minimum in the absence of a recoil kick. 



|Baker et al.] p008] l, & |van Meter et al. 



pOlOt have ail fit distri- 
butions to the data from full NR simulations of BH mergers, with 
results in good agreement. For major mergers (0.25 < q < 1, the 
range we consider in our simulations), high spins (ai = 02 = 0.9), 
and randomly-distributed spin orientations, | van Meter et aL]j2010^ 
find that 69.8% of recoil kicks will be above 500 km s" 1 and 35.3% 
will be above 1000 km s" 1 . If the BH spin magnitudes are also ran- 
domly distributed (0 < 01, a,i < 1), the these fractions are 41.6% 
and 13.5%. Our fiducial-mass merger simulations (described be- 
low) have « csc = 689 — 1248 km s" 1 , and we also use low-mass 
simulations with V csc 3.S low as 338 km s 1 . The NR results indi- 
cate that there is a substantial probability of BHs in these models 
receiving kicks up to ~ v csc . In our high-mass simulations, how- 
ever, v CBC = 1737 — 2555 km s _1 , so the probability of recoil 
events with ~ v asc is much lower in these models. 

We also slightly modify the treatment of BH mergers in the 
code. In the standard GADGET prescription, the BHs merge when 
they are separated by less than a gravitational softening length 
(a scp < -Rsoft) and have a relative velocity v IC \ < 0.5 c s , where 
c s is the local sound speed. If, as can happen in gas-rich merg- 
ers, the central escape velocity is increasing rapidly when the BHs 
merge, small variations in t mrg can result in very different BH tra- 
jectories (see Fig. [5J. In order to reliably compare results from 
different simulations, we therefore force the BHs to merge at a 
given time that is predetermined as follows. We define the "coa- 
lescence time" i coa i as the earliest time after which the BH binary 
is tight enough that it could plausibly merge. Given the numeri- 
cal uncertainties near the spatial resolution limit, we generously 
choose t coa i to be the time when the BH separation falls below 
ciscp = 10 -R so ft (and v TC i < 0.5 c s , as before). For reference, 



-Rsoft 



> pc in our fiducial-mass simulations (see below). Then, 



restarting the simulation slightly before t coa i, we force the BHs 
to merge at a predetermined time t mrg > t coa i. t mTg can then be 
varied as a free parameter, allowing us to systematically probe the 
range of possible BH merger times. In practice, this is only neces- 
sary for nearly equal-mass, gas-rich mergers in which v csc varies 
significantly throughout the merger. In other cases the results are 
insensitive to the choice of merger time, so £ mrg is simply chosen 
to be t coa j. 

For the majority of our simulations, we use the same to- 
tal mass, 1.36 x 10 12 Mq, for the primary galaxy, and scale 
the secondary to yield the desired mass ratio. The exception is a 
small subset of simulations with lower and higher total mass, de- 
scribed below. For all of our "fiducial mass" simulations, we use a 
gravitational softening length of Rsoft = 80 pc for baryons and 
Rsoft, dm = 240 pc for the dark matter (DM). In each galaxy, 
4.1% of the total mass is in a baryonic disk component. The pri- 
mary galaxy has JVhaio = 4.8 x 10 5 DM halo particles, and for 
/gas = 0.4 it has iVdisk = 3.2 x 10 5 disk particles, with equal num- 
bers of gas and star particles. In all other models, iVhalo an d iVdisk 
are set such that the same particle mass resolution is maintained for 
each component (m s tar = 4.2 x 10 5 M , m gas = 2.8 x 10 5 M , 
& m ha io = 5.4 x 10 6 Mq). 

Each galaxy is given a single BH particle, which as previously 
stated, consists of a "seed" BH along with a gas reservoir from 
which an accretion rate is calculated smoothly. We use a small ini- 
tial seed mass for the BH, 1.43 x 10° Mq, in accordance with the 
Mbh - Mbuige ( |Magorrian et al.|1 998) relation, because our initial 
galaxies are pure disks with no bulge component. We find that our 
results are insensitive to the choice of seed mass. The total mass 
of this hybrid particle is the "dynamical" mass of the BH; this is 
the mass that is used for gravitational force calculations. This dy- 
namical BH mass is set to 10 -5 of the total galaxy mass initially. 
This value is chosen such that the seed mass and dynamical mass 
have similar values by the time of the recoil kick, to avoid a dis- 
parity between the mass used for gravitational forces and that used 
for accretion physics of the recoiling BH. To help ensure that the 
BH remains in the center of the galaxy prior to the BH merger and 
recoil kick, we set the accretion drag force equal to its value for 
Eddington-limited accretion prior to the BH merger (though the 
actual accretion rate is still calculated self-consistently based on 
the Bondi-Hoyle formula). Once the BHs merge, the accretion drag 
force is calculated using the actual Bondi-Hoyle accretion rate for 
the remainder of the simulation. 

Our set of fiducial-mass simulations covers a wide range of 
galaxy mass ratios and gas fractions, but in order to probe the ef- 
fects of GW recoil on black hole - host galaxy relationships, as 
we do in § |5.4| we require a larger range in total galaxy mass. We 
have run a sample of seven high and low mass simulations with 
a primary galaxy mass ranging from 0.05 - 20 times the fiducial 
primary mass, which we refer to as Mq. For the lower-mass runs, 
we maintain the same number of particles as in the fiducial simula- 
tions, such that higher mass resolution is achieved. Accordingly, we 
reduce the softening lengths in these simulations by (M/Mo) 1 ^ 3 . 
For the higher-mass runs (Mo x 10 & Mo x 20), an increase in par- 
ticle number by this factor would be prohibitively computationally 
expensive, so we instead increase the particle number by a factor 
of 5 in each case, compromising a factor of 2 and 4 in mass res- 
olution, respectively. In the latter case, we maintain a reasonable 
-AfaH/Mparticic ratio by increasing the initial dynamical mass by a 
factor of 10. 

In order to evaluate the sensitivity of our results to numerical 
artifacts such as the choice of mass resolution, softening length, 
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Figure 1. Fiducial example of a galaxy merger with BHs but no recoil kick 
(model q0.5fg0.3a). The gas distribution at time of BH merger (1.97 Gyr) is 
shown in both the x-y (top) and x-z (bottom) projections, with the density 
scale shown on the right-hand side. The black solid and green dotted lines 
show the pre-merger path of the BHs corresponding to the larger and smaller 
galaxies, respectively. The filled dots denote the BH positions at 500 Myr 
intervals. The black-and-white dot indicates the position of the merged BH 
at the moment of merger. 



and integration accuracy, we have run a number of additional sim- 
ulations in which we systematically vary these parameters. In ad- 
dition to these tests, two of the models included in our results have 
higher mass resolution (and correspondingly higher spatial reso- 
lution) than our fiducial runs: our low-mass merger models have 
the same number of particles as our fiducial-mass models and thus 
have 10 and 20 x higher mass resolution. Higher mass and spatial 
resolution reduce the noisiness of the potential, and varying the res- 
olution changes the detailed structure of the gas, especially in the 
highest density regions. Higher integration accuracy also reduces 
numerical noise. Consequently, the exact trajectory, and hence the 
settling time, of a given recoil event will be affected to some extent 
by these numerical factors. However, the variation in recoil trajec- 
tories for different choices of resolution and integration accuracy 
are small compared to the differences for varying kick speed and 
merger remnant properties. We also define the settling radius of the 
BH to be 4 R ao a to avoid sensitivity of our conclusions to the BH 
motion in the innermost galaxy region where the softened potential 
may dominate. Most importantly, for all choices of mass resolu- 



- 1.000 

2 ° 0.100 

i 0.010 
on 

-g" 0.001 



,— , 10 c 
a 





12 3 4 

Time [Gyr] 

Figure 2. For the same simulation shown in Fig.^(model qO.SfgO.3a), the 
evolution of following quantities throughout the simulation is shown, from 
top to bottom: BH accretion rate (ATbh), global star formation rate, nor- 
malized galaxy gas fraction (/gas//gas,o), BH mass (Mbh), and BH sep- 
aration prior to merger. Vertical lines indicate the time of BH merger. 



tion, spatial resolution, and integration accuracy parameters tested, 
the qualitative behavior of the recoiling BHs and the relative differ- 
ences between simulations are robust. Therefore our main results 
do not depend on these numerical factors. 

Table[T]lists the parameters for the galaxy merger models we 
have constructed. We have tested a total of 62 different galaxy 
merger models in which we vary the galaxy mass ratio (q), the gas 
fraction (/ ga s), the orbital configuration. 55 of these models use the 
"fiducial" primary galaxy mass, and for easier reference we assign 
each of these a name given by q[value]fg[value][orb], where "q" 
denotes the galaxy mass ratio, "fg" denotes the initial gas fraction, 
and each letter orb is identified with a specific orbital configuration. 
The remaining seven models have higher or lower total mass and 
are referred to by q[value]fg[vaIue]M\factor][orb], where "M" de- 
notes the total galaxy mass relative to the equivalent fiducial-mass 
galaxy. 

In all cases, the orbital configuration is specified by six pa- 
rameters, and the galaxies are initially set on parabolic orbits. The 
choice of impact parameter for the first pericentric passage, 7? pe ri, 
and initial galaxy separation, ca, determine the initial orbital en- 
ergy and angular momentum. In our fiducial-mass models, we set 
a; = 143 kpc and R par i — 7.1 kpc, except for orbit "b", which has 
-Rperi = 14.3 kpc. In our high- and low-mass models, a, and R pcl i 
are scaled such that o; is the same fraction (0.625) of the virial ra- 
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dius, i?2oo> and aJRpeii = 20. We denote these orbits with letters 
w - z to differentiate them from the orbits a - k of our fiducial-mass 
mergers. Note that the orbital parameters do not remain constant as 
the merger progresses, owing to energy losses via dynamical fric- 
tion. The angles (81, 4>\) and (82, 4>2) determine the initial orbital 
phase and inclination, respectively, of each galaxy. These angles, 
as well as R pcr i for each galaxy model are given in Table [TJ Our 
fiducial merger orbit (orbit "a") is a "generic" orbit, in that no sym- 
metries exist in the initial orientation angles of the disks. Note that 
many of these orbits are identical to those used in Cox et al. ( 2006 1 
and Robertson et al. (2 006a) . 

For each model, we run a merger simulation in which the BHs 
are not allowed to merge. From this we can determine t coa i, which 
is used to set t mTg . Then, restarting slightly before i mrg , we simu- 
late both a merger with no recoil kick and a merger with Wk/iw = 
0.9. For the models shown in boldface in Table [TJ we also simu- 
late recoil kicks with v k /v csc = 0.4, 0.5, 0.6, 0.7, 0.8, 1.0, 1.1, and 
1.2. Each simulation is run for 2.9 Gyr after i mrg . Finally, using the 
qO.5fgO.3a model, we simulate a small sample of recoil kicks with 
varying kick orientation (#k, 0k), the results of which are discussed 
in §1431 



more compact remnants ^Robertson et al" 2006b| |Dekel & Cox| 



3 MERGERS WITHOUT GW RECOIL 

We begin by examining galaxy merger simulations with no GW re- 
coil kick given to the BHs. Some basic characteristics of a galaxy 
merger simulation with q = 0.5 and / gas = 0.3 (our "fiducial" 
merger, model q0.5fg0.3a) are illustrated in Figs. [TJ & [2] Fig. [TJ 
shows the gas distribution at the time of BH merger, as well as 
the paths of the two BHs prior to merger. The gas distribution 
is very irregular, with distinct tidal streams. The inner region of 
the merger remnant is also lumpy and irregular. In other words, 
the initial disk structure of the progenitors has been destroyed by 
the merger, and the remnant is still highly disturbed at the time 
of BH coalescence. Fig. [2] shows the BH accretion and star for- 
mation throughout the simulation as well as the evolution of the 
BH separation prior to merger. This example illustrates some char- 
acteristics that are generic to our galaxy merger simulations. Si- 
multaneous bursts of star formation and BH accretion occur after 
the first close pericentric passage and at final coalescence (Mihos 
& Hernqui st|1994b| |1996|> as the result of gas inflows caused by 
gravitational torques (Barnes & Hernquist 1991 1996). As we shall 
demonstrate, these gas inflows can greatly influence the dynamics 
of ejected black holes by deepening the central potential and by 
providing additional drag force. Note that because our simulations 
include star formation throughout, the gas fraction at the time of 
BH merger is significantly lower than its initial value; in the exam- 
ple shown, 60% of the initial gas has been consumed by t mrg . By 
the end of the simulation, 2.9 Gyr after the merger, almost 80% of 
the gas has been depleted by star formation, and the BH accretion 
rate is very low. In our other merger models, 47-80% of the initial 
gas is depleted by the time of BH merger, and 62-88% is depleted 
by the end of the simulation. The final black hole mass in this ex- 
ample (2.2 x 10 7 Mq) is also typical of our fiducial-mass mergers. 

Fig. [3] shows gas and stellar density distributions of this fidu- 
cial (qO.5fgO.3a) model and five other merger models. Each plot 
shows gas and stellar density from three different projections. 
These examples are chosen to illustrate the generic nature of the 
morphological features mentioned above. In all cases, the merger 
remnants are visibly disturbed and lumpy. Tidal tails are ubiqui- 
tous, though their size and density varies. Higher- / gas mergers have 



|2006| [Hopkins et al. 2008a Hopkins & Hernqui st|2010a| >. LoweT 
q mergers are less strongly disrupted; a disk-like structure can be 
seen in the x — z orientation of model q0.25fg0.3a. A special case 
is shown in the lower-right corner of Fig. [3] Model qO.5fgO.3h has 
a coplanar orbit with both galaxies rotating prograde to the orbital 
motion. Due to the highly aligned orbit, the disk structure of the 
progenitor galaxies is preserved in this merger remnant; tidal fea- 
tures are visible only from the face-on (x — y) projection. Our other 
"highly aligned" orbits (i, j, & k) result in similar remnant mor- 
phologies; they are included in our suite of merger models mainly 
for comparison, as mergers generally will not have such careful 
alignment. 

As mentioned in § [2] in nearly equal-mass, gas rich merg- 
ers, the potential well in which the BH sits may deepen rapidly 
during final coalescence. Once the merger is complete, the rem- 
nant central potential will become shallower as the central stellar 
region begins to relax. The phase of rapid v csc fluctuation coin- 
cides with the time of BH coalescence and may therefore affect 
recoiling BH dynamics. Fig. [4] shows the evolution of the BH es- 
cape speed (i>csc) for four merger simulations without recoil kicks. 
The solid lines denote t mrg and iw(£mr g ), while the dotted lines 
denote the maximum (post-merger) value of v csc and the time at 
which it occurs. The difference between iw(imr g ) and d CEC)I11 „ 
is small in all except the top left panel, which shows the qlfg0.4a 
merger. It is clear that the escape speed, and hence the trajectory, 
of a kicked BH in this model will depend on whether the kick oc- 
curs at t coa i or at some slightly later time. However, the sharp in- 
crease of u csc occurs only in gas-rich, nearly-equal-mass mergers. 
The other three panels in Fig.|4]show v csc vs. t for a sample of sim- 
ulations with lower q and / gas ; in these models, « csc is much less 
volatile during and after the BH merger. Moreover, the example 
shown in Fig. |4]is not the most gas-rich of our merger models but 
is chosen because of the large difference between v csc (t cos ,\) and 
«csc.m ax , 34%. Larger differences between these quantities yield 
larger uncertainty in our results for BH dynamics when the recoil 
kick is assigned at i mrg = t CO aX- The final four columns in Table [TJ 
give i mrg , Vesc(imr g ), ticscmix, and the fractional difference be- 
tween v csc (t mlg ) and v, 



csc,m a x 



for each merger model. All merg- 
ers with / gas < 0.3 have < 10% difference between « csc (t mr g) 
and i> csc ,max, and mergers with q < 0.5 have < 1% difference. 
Therefore, while our results for nearly equal-mass, gas-rich merg- 
ers are subject to the assumption that the BH merger occurs on a 
short timescale, the results for our other models are insensitive to 
the merger time. We explore further the case of gas-rich, q ~ 1 
mergers in § |4.2| 



4 RECOILING BLACK HOLE DYNAMICS 

4.1 Characterization of Recoil Trajectories 

A universal feature of our recoil trajectories is that they have 
low-angular-momentum orbits, which occurs because the centrally- 
concentrated baryonic component of the galaxy dominates the BH 
trajectories even when they extend far into the halo. Consequently, 
we refer throughout the paper to the trajectories as "oscillations" of 
the BH about the galactic center. We see this clearly in Fig. [5] which 
shows the trajectories of BHs kicked with Uk/tw = 0.4 — 1.0 for 
eight different merger models. Also readily apparent in these plots 
is the variation in BH oscillation amplitude and duration between 
models. For a fixed value of Vk/v Bac , BHs travel further from the 



© 0000 RAS, MNRAS 000, 000-000 



8 Blecha et al. 

q0.5fg0.3a 



I 1 i 1 i Ft 1 I 1 I 




XZ 

• 


— i — 1 1. 1 1 

XZ 

i ■, i ' "r 1 ,.n ' i ' 


■ 1 ' 1 ' 1 ■ 1 ' 1 ' 


y z 


4^ 





qO.5fgO.04a 




q1fg0.3c 




qO.25fgO.3a 



q0.5fg0.3h 



xy 

t 

' I i I i I ' I i I i 


, , ,[;■;. pySfsg^ .1. ] 1 


XZ 


XZ ' 

t?Vy-' ■ 


' 1 ' 1 ' 1 ■ 1 ' 1 ' 
yz 

<* 


. > 2 




Figure 3. Projected gas (yellow/red) and stellar (blue) density distributions for six different merger models, shown at the time of BH merger in each case. The 
black-and-white dot indicates the position of the BH in each panel. The spatial scale of all panels is 86 kpc, and the projection (x-y, x-z, or y-z) is labeled on 
each panel. The models shown are (clockwise from top left, also labeled on each plot) q0.5fg().3a, qlfgf).3a, qlfg0.3c, q0.5fg0.04a, q0.25fg0.3a, & q().5fg0.3h. 



galactic center in mergers with lower q and / ga s- In the qlfg0.6a 
model, BHs kicked with ?; k = 0.7 Vesc (= 878 km s" 1 ) travel 
< 2 kpc from the galactic center, while in the collisionless qlfgOa 
model, the same v^lv csc allows the BH to travel ten times further. 
Note that scaling the kick speeds to tr csc is important for deter- 
mining the trends between galaxy models; iwCimrg) for all of our 
models ranges from 338 — 2555 km s _1 (689 — 1258 km s _1 for 
the fiducial-mass models). The trend toward smaller recoil oscilla- 
tion amplitude for higher / gas occurs because higher gas fractions 
result in more compact remnants that have higher central densities, 



and these remnants also have larger available supplies of gas. Both 
factors contribute to steeper central potentials and increased gas 
drag and dynamical friction. We can also see from Fig.|5]that recoil 
oscillations are slightly larger for lower q, though the trend is more 
evident in higher- / gas remnants. Higher mass ratio mergers create 
stronger perturbations that drive gas more efficiently to the galax- 
ies' central regions, further contributing to the steep potentials and 
high densities. 

We can see these trends more clearly in Fig. [6] The top 
plot shows the maximum galactocentric distance of each orbit in 
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Figure 4. BH v cac vs. time for a sample of four galaxy merger simulations with no recoil kick imparted to the BH. (v csc = y — 2<1?(xbh)-) I n eacn panel, the 
two curves represent v eBC of each BH prior to merger. After the merger, the escape speed of the remnant BH is shown. The solid vertical and horizontal lines 
mark the time of BH merger, t mr g, and the BH escape speed at that time, Dcsc(imrg). The dotted vertical and horizontal lines denote the (post-merger) time of 
the maximum escape speed and its value. Top left panel: qlfg0.4a; top right panel: q().5fg0.4d; lower left panel: q0.333fg0.4a; lower right panel: qlfg0.04a. 



Fig. [5] normalized to the half-mass effective radius of the galaxy 
(i? max /flcft). The bottom plot shows the settling time, t sct tie, de- 
fined as the time when the apocentric distance of the BH oscil- 
lations, relative to the stellar center of mass, falls below 4i? so f t . 
This definition is chosen to avoid following BH trajectories below a 
scale where their oscillations may be dominated by numerical noise 
or gravitational softening. We find that our conclusions are not sen- 
sitive to the exact definition of i se ttie- Note that for some high- 
velocity recoils (for all recoils in the collisionless simulations), we 
have only lower limits on i sc ttic, because in these cases the BH did 
not settle by the end of the simulation. 

In the top panel of Fig. [6] the curves for a given / gas lie nearly 
on top of each other. This supports our assertion that the initial 
gas fraction is more important than the mass ratio in determining 
JJmax/iicff for a given i>k/tw- By normalizing these quantities, we 
essentially isolate the effects of gas drag or potential shape from 
variation in galaxy size or central potential depth. One can easily 
see from Fig. [6] that, especially for lower v^/v eS c, a gas fraction 
> 30% can reduce the amplitude of BH oscillations by up to an 
order of magnitude with respect to a purely collisionless system. 

Note that i? ma xAReff is at least as large for / gas = 0.6 as it is 
for / gas = 0.3, opposite the trend we expect if gas drag dominates 
the BH trajectories. However, galaxies with higher / gas generally 
have higher star formation rates, so systems with initial gas frac- 
tions of 0.3 — 0.6 will have similar amounts of gas by the time 
of the BH merger. Furthermore, as previously stated, higher gas 



fractions yield more concentrated merger remnants (smaller R e g). 
Combined, these effects result in similar R mayi IR c for recoils in 
mergers with initial / gas = 0.3 — 0.6. 

As discussed in § [2] the recoil trajectories may be sensitive 
to numerical parameters such as resolution or integration accuracy. 
In particular, because R e g is of the order of a few kpc, recoiling 
BHs with -Rmax/^cff ^5 0.1 spend substantial time at or below the 
spatial resolution limit. i? max and t BBt tie are therefore approximate 
in these cases. The BH settling times are also susceptible to possi- 
ble inaccuracies in dynamical friction forces due to our finite mass 
resolution, which is likely to result in underprediction of the true 
amount of dynamical friction. However, we can predict that the 
effect of this uncertainty on our results is likely to be small; in 
our lower-mass galaxy mergers, which have 10 x and 20 x better 
mass resolution, BHs kicked near the escape speed still experi- 
ence very little dynamical friction and have very long wandering 
times. Even if dynamical friction is still underestimated in these 
higher-resolution runs, it is clear that in general, recoils near the 
escape speed will result in wandering times of at least a few Gyr, 
and that in gaseous mergers, recoiling BHs with < 0.5 v esc re- 
main within the central kpc of the galaxy. Finally, because we use 
the same mass resolution in all of our simulations (excepting the 
low/high mass runs), the relative trends that we find between dif- 
ferent galaxy models are robust. 

We expect some variation in i? m ax for different merger or- 
bits as well. In the set of eleven different orbital configurations we 
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Figure 5. In each plot, BH recoil oscillation amplitudes are shown for varying kick velocities within a single model. The a;-axis is the time after the BH merger, 
t — tmrg- The color-coded numbers on each plot indicate v^/veac for each curve. Galaxy models shown are, from left to right and top to bottom: qlfg0.6a, 
qO.5fgO.6a, qlfg0.3a, qO.5fgO.3a, qlfgO.la, q0.5fg0.1a, qlfgOa, q0.5fg0a. 
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Figure 6. -R max / Rett (top window) and t se ttle (bottom window) are plot- 
ted versus t>k/i>csc for the eight models whose trajectories are shown in 
Fig. [5] t S ettle is defined as the time after the BH merger at which the apoc- 
entric distance of the BH orbits drops below 4 -R so ft . Arrows denote lower 
limits on t se ttle m cases when the BH has not settled by the end of the 
simulation (defined as t mvg + 2.9 Gyr). 



simulated with q — 0.5 & / g , 



0.3, the central escape speeds 



at i mrg vary from 805 — 1046 km s~ . Since the total mass in all 
of these simulations is constant, the variation in v csc reflects dif- 
ferences in the central concentration, and hence the steepness of 
the potential well. We simulate recoil trajectories in each of these 
with Uk/^esc = 0.9, and find that i? ma xAR c fr = 4 — 20. Based on 
this example, we see that the merger orbit alone can influence the 
remnant central potential well enough to change the amplitude of 
recoil trajectories by a factor of ~ 5. 

Fig. [7] shows two examples of recoil trajectories in three dif- 
ferent projections as well as R(t). The top plot shows a recoil event 
with Uk/^esc = 0.6, which in this model (q0.5fg0.3a) corresponds 
to «k ~ 600 km s _1 . Note that the orbit is highly non-planar; it 
looks similar in all three projections. The orbit is also centrophilic, 
which as we have noted is common to all of our recoil simulations. 
The average orbital period is short, ~ 7 x 10 5 yr, and i sc ttic ~ 400 
Myr, such that the BH settles to the center well before the end of 
the simulation. 

Recoil events with velocities near « CS c produce long-lived (> 
Gyr), large-amplitude BH oscillations. Many of these oscillations 
still have large amplitudes at the end of the simulation, 2.9 Gyr af- 
ter the BH merger, so we have only lower limits on i 8e ttle- The 
bottom plot in Fig. [7] illustrates this type of recoil event, again 
for the qO.5fgO.3a model. The trajectory shown has a kick speed 
of Vk/fW = 0.9, so Uk ~ 900 km s -1 . The amplitude of these 





-0.6 -0.3 0.0 0.3 0.6 
x [kpc] 



-0.6 -0.3 0.0 0.3 0.6 
x [kpc] 




-0.6 -0.3 0.0 0.3 0.6 

y [kpc] 



0.6 

0.4 

0.2 
0.0 



III ' 

ii Viiiiiii 




l|| 




Itfthi 



0.0 0.1 0.2 0.3 



Vv esc = 0.9 



40 








40 


20 








Q. 


20 









Z[k 





20 








-20 


40 








-40 




-40 -20 20 40 
x [kpc] 



-40 -20 20 40 
x [kpc] 




40 


A 


A 


A 


A 


A; 


30 

Q. 












re 20 












10 
























-40 -20 20 40 
y [kpc] 



0.0 0.5 1.0 1.5 2.0 2.5 
t-UfGyr] 



Figure 7. Top four windows: trajectory for the v^/v eBC = 0.6 recoil sim- 
ulation with the q().5fg0.3a model. Bottom four windows: trajectory for 
the Dk/fesc = 0.9 recoil simulation with the q0.5fg0.3a model. In each 
case, the four windows show the trajectories in the x — y, x — z, and y — z 
projections, as well as the BH's galactocentric distance versus time. Note 
the greatly different spatial scales of the two simulations. 



oscillations is ~ 50 kpc, almost two orders of magnitude greater 
than in the lower-^k example. Although the BH only completes 
five orbits before the end of the simulation, we can see that again 
there is no evidence for a preferred orbital plane; the orbit is three- 
dimensional. Because the oscillation amplitude is so large and 
shows no sign of damping out by the end of the simulation, it is 
likely that such an orbit will not return to the galaxy center within a 
Hubble time, and certainly not within the mean time between major 
galaxy mergers. Such BHs essentially can be considered "lost" to 
the galaxy. If in other cases the BHs do return to the galactic center 
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after several Gyr and another galaxy merger has taken place in the 
meantime, the returning and new BHs could form another binary. 

In recoil events with v > v csc , the rapidly recoiling BH can 
only take a small amount of gas with it; only gas that is orbiting 
the BH with v OI h > «k stays bound to the BH when it is ejected. 
The accretion timescale and luminosity of these ejected disks is 
discussed in § |5.2| Once this supply is exhausted, the naked BH is 
unlikely to become active again, as this would require a serendip- 
itous passage through a dense gas cloud, as well as a sufficiently 
low relative velocity that substantial accretion could occur. (Recall 
that Bondi-Hoyle accretion scales as «gg.) These escaping BHs 
are therefore likely to wander indefinitely, undetected, through in- 
tergalactic space. 

4.2 BH Merger Time 

As illustrated in Fig. [4] i> csc can increase significantly during coa- 
lescence in equal-mass, gas-rich mergers. To better understand how 
this may affect BH trajectories, we have run a series of simulations 
with the BH merger (and recoil kick) occurring at sequentially later 
times using model qlfg0.4a. The v asc evolution of model qlfg0.4a 
is replotted in the top window of Fig. [8] with the fiducial merger 
time, t mrg = t coa i, marked by the solid black line and the subse- 
quent merger times tested marked by black dashed lines. We have 
tested simulations with each of these merger times twice, once 
with iik = 1000 km s _1 and once with Wk/fesc^mrg) = 0.9. By 
conducting the second experiment we are normalizing to the evolv- 
ing depth of the central potential well, which helps us to separate 
the gravitational and drag effects on the trajectories. 

Fig.[5]shows Jlmax for the BH orbits in each of these varying- 
t mrg simulations, with the = 1000 km s simulations in the 
middle panel and the fk/fesc(tmr g ) = 0.9 simulations in the bot- 
tom panel. As v CBC increases from 1100 km s _1 at i coa i to al- 
most 1500 km s _1 at its maximum, the amplitudes of BH trajec- 
tories from constant-Uk recoil events decrease by more than two 
orders of magnitude. However, the bottom panel of Fig. [8] shows 
that when is set to a constant fraction of the escape speed, 
the BHs kicked from deeper potential wells (larger v esc ) still have 
smaller-amplitude trajectories, but the variation is much less se- 
vere than when «k = constant. The variation that persists even for 
fixed Vk/v esc can be attributed to the enhanced gas drag and dy- 
namical friction, as well as the steeper central potential well, that 
result from the rapid formation of the central dense cusp. At very 
late merger times, 100 — 500 Myr after t C oal> the opposite trend 
is seen. The amplitude of constant-?^ trajectories increases by al- 
most an order of magnitude from the minimum at i coa i + 70 Myr 
to the value at t coa _i+ 500 Myr. This is partly because the central 
potential becomes shallower as the merger remnant begins to relax. 
The same trend occurs at late times in the constant-Uk/« C sc simula- 
tions, because the potential well also becomes less steep, and gas 
drag becomes less efficient as more of the gas is consumed in star 
formation. We note that merger times up to £ coa i + 500 Myr are 
included here only for illustrative purposes; in reality, such delayed 
mergers are probably unrealistic in gas-rich systems. 

If the BHs are able to merge rapidly, before the central po- 
tential reaches its maximum depth, another interesting effect can 
occur. During the period of rapid central potential evolution, the 
timescale on which v csc increases is typically much shorter than 
the timescale for a recoiling BH orbit to decay via gas drag or dy- 
namical friction. Consequently, on its first several pericentric pas- 
sages through the central region, the BH encounters an increas- 
ingly deep potential well, and its velocity is actually boosted by 
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Figure 8. Top window: v eBC is plotted as a function of t— t CO al for the 
qlfg0.4a merger simulation. The portion shown here is a zoom-in of the 
top-left window in Fig. [5] The solid vertical line denotes t coa l, and the 
dashed vertical lines denote the subsequent merger times used to test vari- 
ation in t mrg . Middle window: -Rmax/^eff is plotted versus BH merger 
time for the set of delayed-merger simulations done using model qlfg0.4a. 
Each recoiling BH in these simulations is assigned a recoil kick of 1000 km 
s — 1 . Bottom window: Same as middle window, but for a set of simula- 
tions in which the kick speed in km s _1 varies, but is held to a con- 
stant value of v^/v eBC ; in every case v^/Venc = 0.9. The x-axis ranges 
from i mrg = tcoal to t mrg = t coal + 500 Myr. 



a small amount, such that the ratio VBn/v eac is roughly constant on 
subsequent pericentric passages. Because only vbh in creases, not 
vbh/v csc , the BH's galactocentric distance is not boosted beyond 
its initial maximum, i? max , and once the rapid potential evolution 
ceases, drag forces become dominant and the velocity begins to de- 
cay. This effect is by design short-lived, and it is strongest when 
the BH orbital timescale is short, such that the BH orbit is domi- 
nated by the dynamics of the inner galaxy region. These velocity 
boosts are therefore mostly an interesting illustration of the sensi- 
tivity of recoil dynamics to the BH merger time, although we will 
see some implications of this effect for the recoiling AGN lifetimes 
calculated in § |5.3| 

The primary implication of our BH merger time analysis is 
that in nearly equal-mass, gas-rich mergers such as the example 
shown in Fig. [8] a delay of even a few x 10 7 yr in the BH merger 
may significantly reduce the amplitude and duration of recoiling 
BH oscillations. Based on results from numerous hydrodynamic 
simulations of binary BH inspiral (e.g., |Escala et al.||2005] |Dotti| 
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Figure 9. Top window: distribution of -R max for our set of (6, (p) 
5^ (0,0) simulations, normalized to i? max (0, 0) for each correspond- 
ing Ufc/uesc value. The set of simulations includes runs with tik/i'csc = 
0.6, 0.8, & 0.9. The black solid line is the distribution of simulations with 
the q0.5fg0.3a model, and the red dashed line corresponds to the qlfg0.3a 
model. Bottom window: the distributions of the average orbital period 
are shown for the same simulations as in the top window, normalized to 
(P orb (0,0)>. 

|et al.|20 07), we have good reason to believe that BHs in gas-rich 
systems merge efficiently in general, but the merger time of a given 
BH binary is impossible to predict with this precision. By assum- 
ing i mrg = t coa i in our simulations, we exclude the possibility of 
delayed mergers, but we also take steps to ensure that this has a 
minimal effect on our results. First, only a small number of our 
simulations have large variation in v csc . In these cases, we see by 
comparing Figs. [6] & [8] that the variation in oscillation amplitudes 
for a given Vk/Uesc due to evolving central density (less than a factor 
of 10) is generally much less than the variation for different values 
of Wk/fcsc (up to three orders of magnitude). All of the analysis in 
this study is presented in terms of v^/v csc . Therefore, the dramatic 
suppression of recoil trajectories with fixed «k, while an important 
result in itself, does not affect our other conclusions. 

4.3 Recoil Kick Direction 

Most of our recoiling BHs have > 500 km s _1 , in which case 
their kicks will be oriented out of the orbital plane of the binary 
( |Campanelli et al.|2007a|b)|Lousto & Zlochower|2009| l. However, 
because the BH binary inspiral occurs on sub-resolution scales in 
our simulations, the binary orbital angular momentum vector can- 
not be predicted. Furthermore, the galactic disk structure of the 
merger progenitors is disrupted by major mergers (see Fig. [3}, 
except when the galaxies have nearly coplanar orbits. In generic 
mergers, there is no obvious reference direction in which to ori- 
ent the recoil kick with respect to the progenitors' initial structure. 
Clearly, on certain recoil trajectories (i.e., along tidal tails or central 
clumps) a BH will encounter denser gaseous regions than on oth- 
ers. However, it is unclear a priori whether these density variations 



are sufficient to cause significant variation in BH recoil trajectories. 
In most of our simulations, we simply orient the recoil kick along 
the z-axis of the simulation, (#k,<^k) = (0,0), in the coordinate 
system with respect to which the initial galaxy orbits are assigned. 
Here, we examine a subset of simulations with varying kick orien- 
tation to determine the sensitivity of our results to this choice. 

We test a sample of different kick orientations with three 
different values of v^/v csc (0.6, 0.8, & 0.9) in our qlfg0.3a and 
q0.5fg0.3a merger models. As a gauge of how the BH trajecto- 
ries differ among this sample, we calculate the maximum galac- 
tocentric distance (i? max ) and average orbital period ((P rb)) of 
the BH in each simulation. Fig [9] shows the distributions of these 
quantities, normalized to _R max and (P or b) of the (9^, <j>k) = (0, 0) 
simulations for each corresponding v^/v csc value. The qO.5fgO.3a 
distributions peak slightly below 1.0 and the qlfg0.3a distributions 
peak slightly above 1.0, but the combined histogram is distributed 
fairly evenly around unity. In other words, the (#k,0k) = (0, 0) 
kick orientation is not a "special" direction. The overall variation 
in -RmaxARmax(0, 0) for each model is less than a factor of 4. In 
contrast, Fig. [5] shows that i? max is almost two orders of magni- 
tude higher for Vk/iW — 0.9 than for 0.6. Clearly, the amplitude 
and duration of BH oscillations depends much more strongly on 
the kick magnitude than its direction, and we need not concern our- 
selves too much with our choice of kick orientation. 



5 BLACK HOLE ACCRETION AND FEEDBACK 

5.1 Bondi Accretion and AGN Lifetimes 

Fig. [TO] compares the BH accretion and star formation histories 
for simulations of the q0.5fg0.3a merger model with w k = 
and «k/«csc = 0.5 — 0.8. The top windows in Fig. |10| show that 
even kicks with Vk/v esc = 0.5 — 0.6 produce qualitatively differ- 
ent AGN lightcurves (as inferred from A/bh) than in the zero-kick 
case, although the final BH masses are similar. In these cases, the 
peaks in accretion rate at coalescence (t « 2 Gyr) are slightly lower 
due to the sudden removal of the BHs from the highest-density re- 
gion. However, the BHs remain within 1 kpc of the galactic center, 
and accretion continues through the end of the simulation at higher 
rates than in the zero-kick case. Thus, although the total mass ac- 
creted is about the same for Vk/i>esc = — 0.6, the active lifetime 
is longer when a recoil kick occurs. 

In the Uk/i>esc = 0.7 simulation, the BH accretion cuts off 
more sharply at i mrg and is subsequently highly variable, due to 
the larger-amplitude, longer-lived BH oscillations. Unlike the BHs 
with lower-velocity recoils, this BH attains a noticeably lower final 
mass than its zero-kick counterpart. 

For recoil kicks > 0.8 v csc , the BH accretion cuts off sharply 
at tmrg, and the BH subsequently accretes very little gas as it trav- 
els on large, long-period orbits that extend well into the halo. These 
high- velocity recoil events therefore reduce the AGN lifetime. Ad- 
ditionally, the BH mass is essentially frozen at its value at t mrg , 
which is about half that of the stationary BH mass. If the BH does 
return to the galaxy center after several Gyr, most of the galaxy's 
gas will have been consumed in star formation, so the BH will re- 
main undermassive. In the galaxy merger models we have used, the 
BH mass deficits for high-velocity recoils versus no recoil range 
from a factor of ~ 1 (almost no deficit) to ~ 5. These BH mass 
deficits and their implications are discussed in § |5.4| 

To quantify the amount of enhancement or reduction in AGN 
activity owing to GW recoil, we calculate AGN lifetimes (£agn) 
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Figure 10. The same quantities for model q0.5fg().3a are plotted as in Fig. [2] but here they are shown for simulations with v^/vcsc = 0.5 — 0.8 in addition to 
the = case. The red lines in the Mbh plots (top row) correspond to the accretion rate calculated from our ejected-disk model as described in 
the t];/0 plots, the v^ = curves are drawn again with dashed lines, except for the A/bh plots where the dashed lines are omitted for clarity. 



5.2 



In 



for both stationary and recoiling BHs in six different galaxy mod- 
els for which we have varied v^. Fig. |1 1| shows the results, with 
colored lines indicating £agn for recoiling BHs as a function 
of and black horizontal lines indicating 4agn for a station- 
ary BH in the corresponding merger model. In order to differ- 
entiate between bright AGN and low-luminosity AGN, we calcu- 
late t agn assuming three different definitions of an "active" BH. 
The thick solid lines in the figure are calculated assuming the BH 
is an AGN when Lboi > 10%Z/Edd, and the thin solid lines as- 
sume Lboi > 3% I/Edd- Because these quantities depend on the BH 
mass, we use an absolute luminosity value for our third AGN def- 
inition: Lboi > Lmin (denoted by the dashed lines in Fig.^TJi. We 
choose L m i n = 3.3 x 10 9 L , which is 0.1 Lsdd for a 10 6 M Q BH. 



In these models, the low-luminosity AGN lifetime is gener- 
ally enhanced by moderate-velocity recoil events, and the bright 
AGN lifetime is generally reduced. Although the details depend on 
the model parameters and kick speed, in all models, Jagn falls off 
steeply at high v^/Veac, when the BH is ejected far enough from the 
center on large enough orbits that Bondi accretion becomes inef- 



fective. This corresponds to the accretion behavior seen in the last 
column of Fig.[l0](u k /i; csc = 0.8). 

For smaller vw, such that the BH is confined to the central few 
kpc, the dependence of £agn on the nature of the merger remnant 
is more complex. We examine first the relatively gas-poor models 
(those with initial / gas = 0.1). In the qlfgO.la model, £agn is en- 
hanced only for v^/v CEC < 0.5 and falls precipitously at higher v^. 
In the q0.5fg0. la model, recoil events always reduce the AGN life- 
time. Merger remnants with less gas have smaller central gas reser- 
voirs to fuel the BH during pericentric passages, and their shallower 
central potentials result in larger BH orbits, which have fewer peri- 
centric passages. Consequently, in gas-poor remnants, t agn will 
generally be lower for recoiling BHs than for stationary BHs. 

In gas-rich merger remnants, £agn is more readily enhanced 
by GW recoil. Stationary BHs in gas-rich mergers typically expe- 
rience a bright AGN phase due to the rapid inflow of cold gas to 
the central region, which is terminated when the AGN feedback 
energy heats the surrounding gas enough to quench the BH fuel- 
ing. Fig. [TTJ shows that although GW recoil shortens this bright 
AGN phase, the low-luminosity t agn is enhanced for recoil kicks 
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Figure 11. AGN lifetimes (4agn) of recoiling BHs with v^/v^c = 0.4 — 0.9, compared to <agn f° r = in six different galaxy models. In each plot, 
the black lines indicate the stationary AGN lifetimes and the colored lines indicate the recoiling BH lifetimes at each value of ife/t^sc- The models shown 
are (from top left): qlfg0.6a (red), qlfg0.3a (blue), qlfgO.la (green), q0.5fg0.6a (magenta), q0.5fg0.3a (cyan), and q0.5fg0.1a (gold), as indicated in the plot 
labels. The different line types denote f agn according to three different definitions of an "active" BH. For the thin and thick solid lines, t agn is defined as the 
time for which the BH luminosity Lbol > 3% and 10% Z/Edd> respectively. The dashed lines show *agn f° r which Lbol > ^min, where L m j n = 3.3 X 10 9 
Lq, or O.lLBdd for a 10 6 Mq BH. In all cases Lbol is calculated from the Bondi-Hoyle accretion rate as described in the text. No cuts have been applied 
to *agn based on the BH settling time, t sc ttle> so in some cases *agn includes accretion that occurs after the BH settles back to the galactic center. 



as high as Vk/v esc = 0.7. In these cases, the recoiling BHs are on 
tightly-bound orbits and can accrete gas during their numerous pas- 
sages through the central dense region. Because they remain in mo- 
tion, however, their feedback energy is spread over a much greater 
volume than for stationary BHs. As a result, the gas encountered 
by the recoiling BHs never heats up enough to completely cut off 
the BH fuel supply, and the BH has a longer lifetime as a low- 
luminosity AGN. In contrast to the traditional "feast or famine" 
model of merger-triggered AGN fueling, these short-period recoil- 
ing BHs are more inclined to "nibble". 



We can see evidence of this displaced AGN feedback through 
its effect on star formation. The bottom row of Fig. [10] shows the 
depletion of gas throughout each simulation, normalized to the ini- 
tial gas fraction, and the dotted lines show the «k = curve in each 
of the other plots. It is evident that slightly less gas remains at the 
end of each higher-u k simulation. In the merger with a stationary 
BH, this gas is accreted or expelled from the galactic center via 
feedback, but when the central BH is removed, the gas continues 
be consumed by star formation. Although in this example the star 
formation rate enhancement is small and difficult to discern from 



© 0000 RAS, MNRAS 000, 000-000 



16 Blecha et al. 



examining the SFR evolution alone, the effect can be much larger 
in higher-q, / gas mergers. We explore in § |5 .5| the implications of 
enhanced SFRs for our equal-mass, gas-rich mergers. 

Although we have found that i agn may be enhanced by GW 
recoil for low-luminosity AGN, we emphasize that the bright AGN 
phase (L > 10%1/Edd) is always shorter for recoiling BHs than for 
stationary BHs, or at best roughly equal. This is a direct result of 
the peak accretion episode being disrupted at the time of the recoil 
kick. The important consequence of this is that the total accreted 
BH mass is never enhanced by GW recoil in our simulations. For 
moderate recoil speeds, the BH mass deficit relative to a stationary 
BH is generally slight, but we find that the effect of recoil is always 
to decrease the final BH mass. The implications of this finding are 
discussed in § |5.4| 

Finally, we note a possible observational consequence of 
longer AGN lifetimes. Fig. [TU] clearly shows that when — 0, 
the two main accretion episodes are associated with the periods of 
rapid star formation. This means that the AGN luminosity may be 
difficult to distinguish from the starburst luminosity, and the AGN 
may also be dust-obscured for at least part of its active phase. If 
GW recoil allows the AGN to remain active long after the starburst 
is complete, such an AGN might be more readily observable, be- 
cause it would no longer be competing with the starburst luminosity 
and may also be less dust-obscured. 

5.2 Ejected Gas Disk 

We have seen that moderate-velocity recoils in gas-rich remnants 
may allow efficient accretion from ambient gas. However, a BH 
moving rapidly or in a low-gas-density region is unlikely to sweep 
up much gas from its surroundings, and we would expect the prob- 
ability of observing this type of active BH to be low. Indeed, in our 
simulations, the Bondi accretion rate of the BH rapidly declines 
following a large recoil kick (see Fig. 1 10|>- The possibility not ac- 
counted for in our simulations is that, if the BH is surrounded by 
an accretion disk at the time of recoil, some amount of this gas may 
remain bound to the recoiling BH. We therefore implement an ana- 
lytic model to calculate the BH accretion rate from this ejected gas 
disk, based on the BH mass and accretion rate at the moment of 
recoil. 

The BH carries with it the part of the accretion disk where 
forb ^ Vk', i-e., it will carry a disk with radius approximately 
R c j ~ GMbh/v^ 2 . This is a reasonable approximation if the BH 
is kicked directly out of the orbital plane of the disk, as is the case 
when Wjj > 500 km s _1 and the BH orbital plane is aligned with 
that of the disk. If the accretion rate from this ejected disk is high 
enough, such a system could potentially be seen as an offset AGN 
via either resolved spatial offsets or offset spectral lines. BL08 pre- 
dicted accretion timescales of < 1 — 10 Myr for the type of re- 
coiling BHs we consider here (Mbh ~ 10 6 — 10 s M0; «k ~ 
400 — 1000 km s -1 ), assuming a constant- M model. Here, we im- 
prove on these estimates by calculating a time-dependent accretion 
rate for the ejected disk. 

We assume the innermost region of the gas disk is well- 
modeled by a viscous, Keplerian, thin disk (specifically, an a-disk, 
|Shakura &"S unyaev ( 1973 1). The evolution of the disk surface den- 
sity, E, can be described by a diffusion equation: 



R 



1 2 d 
dR 



vYtR 



1/2 



(1) 



where v is the kinematic viscosity (see, e.g., Lynden-Bell & Pringle 
1 1974} |Pringle|1981| >. The steady-state solution is recovered when 



the left-hand side of the equation vanishes, but this is not the case 
for our ejected disk, which is no longer fed at its outer radius. Exact 
analyti c solutions to Eq. {T]> can be f ound when v scales only with 
radius ( Lynden-Bell & Prin gle|1974| l, but in the a-disk model, v is 
a function of both R and E. Pringle ( 19911 has derived self-similar 
solutions to this equation for viscosity scaling as v oc E m i?". 
These solutions assume a disk with no inner radius, so for our case 
we must neglect the finite inner disk radius at the innermost sta- 
ble circular orbit of the BH (-RiscoX as well as the possibility of 
a circumbinary gap in the disk that persists after the BH merger. 
This should not affect our results, however, as the ejected disk will 
be much larger than the ISCO (Rcj/Risco > 10 3 ), and the gap 
is expected to refill on a timescale of years (Milosavljevic & Phin- 
ney |2005l|Tanaka & Menou|20T0l >. For an a-disk in which viscos- 
ity is assumed to scale with gas pressure and the opacity equals 
the electron-scattering value (k cs « 0.4 cm 2 g _1 ), the viscosity is 
given by 

v = CRY? lz , 

4 / 3 , 



(2) 
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In this equation fca is Boltzmann's constant, a is the thin-disk vis- 
cosity parameter, fj, is the mean molecular weight, m p is the proton 
mass, and <tsb is the Stefan-Boltzmann constant. This corresponds 
f | Pringle|(l99T) with m = 2/3, n = 1 



to the self-similar solution of 



(see also Canniz zo et al.|1990| 
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The arbitrary constants Eq, to, and Rg satisfy 
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Integrating Eq. ([3} yields a time-dependent expression for disk 
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From this we obtain the time-dependent accretion rate: 
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We determine Eq, to, and Rq using Eq. |4| along with the con- 
ditions that Md(t') = M ej and M d (t') = M mrg , where t' = 



to + t s 



- t n 



tsim is the time in the simulation, M e j is the 



ejected disk mass, and Af mrg is the Bondi accretion rate calculated 
in the simulation at the time of merger. The mass evolution of the 
disk is therefore given by: 

-3/16 

(7) 



Md(i') 



Md(t') = M n 



to 
to = 



-19/16 



3 M c 



16 M n 



(8) 



(9) 



We can calculate the initial radius and mass of the ejected disk 
using the multi-component, self-gravitating disk model of BL08. 
The calculation takes as input the kick velocity, BH accretion rate, 
and BH mass at the time of the merger. Note that because we use 
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the accretion rate calculated in the code as input, the normaliza- 
tion of the ejected-disk accretion rate is subject to the same as- 
sumptions as is the Bondi-Hoyle accretion rate (see § |2j. We re- 
fer the reader to BL08 for the details of the self-gravitating disk 
model, but we note that for the BH masses and accretion rates 
considered here, R c j lies in the region of the disk that is self- 
gravitating but still strongly dominated by the BH potential (i.e., 
Zone II of the BL08 model). Quantitatively, R e j for our fiducial- 
mass simulations is between lo 3 - 9-5 - 7 Schwarzschild radii, or 
~ 0.02 — 0.6 pc. The corresponding mass of the ejected disk, M e j, 
is 10 3 ' 5 " 7 ' 5 M (lo 4 6 " 6 3 M for our fiducial-mass models). 

We use these quantities as initial conditions for the ejected 
disk model outlined above. The timescale for the accretion rate (and 
mass) decay, to, ranges between io 5 - 6-7 - 6 yr. to depends on M e j 
and M mrg , which in turn depend on Vk and the galaxy model pa- 
rameters. Specifically, to decreases with increasing «k, / gas , and 
q. After a time to, the accretion rate has dropped by a factor of 
2-19/16 w o 44 ( note mat tne tmle m Eqns.[9 is defined such that 
t! = to at i mrg ). By the end of our simulations, about 55 — 80% of 
the disk mass has been consumed. This corresponds to an increase 
in the BH mass of AM = 0.3 - 6% (0.3 - 3% for v k /v csc = 0.9 
simulations only). As expected, these values of AM are somewhat 
lower than the result of BL08, AM w 10% Mbh. BL08 consid- 
ered kick speeds as low as 100 km s _1 ; their result also assumed 
that the entire disk mass was accreted and that the accretion was 
Eddington-limited at the time of BH merger. 

Fig.[l0]demonstrates that except in high-fk simulations where 
the average orbital timescale is long (vk > 0.8 v csc ), the ejected- 
disk accretion is typically insignificant compared to the Bondi ac- 
cretion of the oscillating BH. Of course, the BH may still carry an 
accretion disk during this phase, but the accretion will be more ac- 
curately characterized by the highly-variable Bondi rate calculated 
from the local density and sound speed than by an isolated-disk 
model. Therefore, our analytic ejected-disk model becomes impor- 
tant only for high-velocity recoils in which the Bondi accretion rate 
drops sharply after the kick. Below, we discuss the AGN lifetimes 
resulting from the combined accretion models, as well as the feasi- 
bility of observing offset AGN. 



5.3 Offset AGN 

In Fig. ^] we show 2 examples of the AR, Av phase space oc- 
cupied by recoiling BHs. The left and right panels show simula- 
tions of the q0.5fg0.3a model with Vk/v csc = 0.6 & 0.9, respec- 
tively. AR is the binned spatial offset of the BH from the stel- 
lar center of mass, and Av is the binned kinematic offset from 
the stellar center-of-mass velocity. The contour shading shows, 
for each AT?, Av bin, the amount of time spent at that location 
(At, top windows), the total AGN energy output (ASboi, middle 
windows), and the average bolometric luminosity (Lboi) (bottom 
windows). Lboi is determined at each timestep from the relation 
Lboi = e r ad Mbh c 2 , where the radiative efficiency e ra( j is 0.1 un- 
less Mbh < O.OlMEdd, in which case the system is consid ered ra- 
diatively inefficient, and e ra d = 0.1 AfeH/(0.01MEdd) iNarayan 
|& McClintock|2008) . 

The Uk/«csc = 0.6 recoiling BH, which settles back to the 
galaxy center in ~ 400 Myr, fills a corner of the phase space. 
The At contours peak at low velocities, as these represent apocen- 
tric passages where the BH spends the most time. (Lboi) is above 
lO 8 I/0 throughout the simulation, with peak values > W 11 Lq. 
Unlike the At contours, (Lboi) peaks at high velocities, because 
the BH has the highest accretion rates immediately after the kick 



and during pericentric passages. The total bolometric energy out- 
put per A_R, Av bin, which is the product of these two distribu- 
tions, accordingly peaks at both high AR, Av and low AR, Av. 
Instead of emitting all of its AGN feedback energy from a single 
central point in the galaxy, the BH distributes this energy through- 
out the central kpc, as discussed in § |5.1| Consequently, this BH 
has a longer AGN lifetime than the v k — BH in the same merger 
model (see Fig.|l 1}. 

In contrast to this example, the Vklvesc — 0.9 recoiling BH oc- 
cupies a narrow track in the phase space, as it completes only a few 
orbits and does not experience much drag. The At and (Lboi) con- 
tours in this example illustrate an inherent difficulty in observing 
rapidly-recoiling BHs: the brightest luminosities occur just after 
the kick, when the BH velocities are highest, but this is also where 
the BH spends the least amount of time. Still, the BH has luminosi- 
ties > 10 8 Lq for AR up to 10 kpc. Further analysis is required to 
determine the observability of the recoiling BHs in both of these 
examples, as well as those in our other simulations. 

Fig.[T3]shows the AGN lifetimes, £agn, in six merger models 
for which we have varied v k . AGN lifetimes in the top two windows 
include only the Bondi accretion rate, while £agn in the bottom 
two windows includes the accretion from our ejected-disk model as 
well. Unlike the lifetimes in Fig.[TT| here we also apply a cut to the 
AGN lifetime such that only activity that occurs while the BH is re- 
coiling is shown: £agn < tscttic- We calculate t agn using each of 
the three AGN definitions described previously. The solid lines in 
the figure correspond to AGN defined as Lboi > (3%, 10%)LEdd, 
and the dashed lines correspond to Lboi > L m i n . The shaded re- 
gions of the plot show the range of values between these definitions. 
As before, we choose L m i n = 3.3 x 10 9 Lq = 1.3 x 10 43 erg s _1 . 

A notable feature in Fig.[T3]is that the peak AGN lifetimes for 
recoiling BHs are quite long, > 100 Myr, even for kick speeds near 
1000 km s _1 in some cases. Note also that in the lower- / gas merg- 
ers, £agn peaks at the lowest Vk/v BSC , but in the gas-rich merg- 
ers the peak value is v k /v csc = 0.7 — 0.8. This reflects the fact 
that we have defined 4agn < tsettie in these plots. In the gas- 
rich mergers, gas drag quickly damps out the low-?;k trajectories, 
and £ agn is limited by tgettie- For higher Vk/v csc , gas drag is less 
efficient and £agn is instead limited by the gas supply available to 
the BH. 

In all cases, the AGN lifetimes that include only Bondi ac- 
cretion (Fig. ^] top windows) fall sharply to < 1 Myr for 
large v k /v csc . This demonstrates that Bondi accretion is inefficient 
when the BH is kicked far from the dense central region. How- 
ever, when we also allow for BH accretion from its ejected disk 
(Fig. [13] bottom windows), the ejected-disk accretion enhances 
the active lifetime for high Vk/v esc recoils by more than an order 
of magnitude. In the examples shown, £agn for low-luminosity 
AGN is 3 - 14 Myr even for kicks as high as 1.2 v csc . For even 
higher kick speeds, £agn will slowly decay as the trend in the 
plot indicates, owing to the decreasing ejected disk mass. The rel- 
atively slow decay of 4agn with increasing i)k occurs partly be- 
cause, as opposed to a constant- M accretion model, here t agn de- 
pends on the fraction of disk accreted before the BH luminos- 
ity falls below the "AGN" limit. In addition, the disk mass in the 
self-gravitating regime has a fairly weak (r 1 ^ 2 ) dependence on ra- 
dius, so M c j oc 1/wk- Recall, however, that the probability of re- 
coil events with kick velocities ^> 1000 km s _1 is quite small 



(|Schnittman & Buonanno||2007| |Campanelli et al.||2007a| |Baker| 
|et al.|2008||Lousto et al.|2010b||van Meter et al.|2010| >; this more 
than the BH's fuel supply limits the chances of observing a recoil 
event with Uk 2> i> OS c. 
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Figure 12. Contour plot of recoiling BH quantities At (top panels), ABbol (middle panels), and (Lbol) (bottom panels) in AR, Av space. (Lbol) is the 
average bolometric BH luminosity per bin. The left panels are calculated from the q0.5fg0.3a merger simulation with v^/veBc = 0.6, and the right panels are 
for the same merger model but with v^/v CBC = 0.9. 



These findings illustrate the importance of BH accretion from 
the ejected disk for the lifetimes of rapidly-recoiling AGN. Also 
crucial, however, is the time-dependent nature of this accretion. If 
we were to assume a simple constant- M model for the ejected-disk 
accretion, we would greatly overpredict the lifetime of the bright 
AGN phase, and in some cases we would underpredict the lifetime 
of the low-luminosity phase. We can calculate the discrepancies 
between £agn for a constant-M model and the lifetimes shown 
in Fig. [13] restricting the comparison to Vk/v CBC > 0.9 to capture 
only the simulations in which ejected-disk accretion is the domi- 
nant mode. MBH(imrg) ~ 10% — 50%MEdd for all the models 
shown in Fig. [13] except q0.5fg0. la. If we assumed the accretion 
continues at this rate for a time Af c j/A/BH(tmrg) in these models, 
we would overpredict £agn for bright (L > 10%Z/Edd) AGN by 
a factor of ~ 4 — 40. Conversely, £agn for L > L m i n would be 
underpredicted by a factor up to ~ 7. In the q0.5fg0. la model, 
MBH(tmrg) ^ 5%MEdd- In this case, £agn for L > 3%Z/Edd 
would be overpredicted by a factor of ~ 10 — 40, depending on 
kick speed, and £agn for L > L m i n would also be overpredicted 
by a factor of ~ 2 — 4. Therefore, it is clear that inclusion of more 



realistic, time-dependent model for ejected-disk accretion, such as 
the one presented here, is important for calculating the lifetimes and 
luminosities of rapidly recoiling AGN. Similarly, for lower-velocity 
recoils, the BH + ejected-disk system clearly cannot be treated as 
isolated, because accretion from ambient gas is the dominant mode 
of AGN fueling. 

To determine whether these recoiling AGN could be distin- 
guished from their stationary counterparts, we need more informa- 
tion than simply the active lifetimes. Recoiling AGN could be iden- 
tified as such if during their active phase they are either kinemati- 
cally or spatially offset from their host galaxies; we will consider 
both possibilities. 

The gas carried along with the recoiling BH is expected to 
include at least part of the BH's broad-line (BL) region, because 
BL clouds reside deep within the potential well of the BH. The 
narrow-line (NL) region lies much further from the BH, so in gen- 
eral this gas will not be bound to the recoiling BH and will remain 
at the redshift of the host galaxy. A kinematically offset AGN could 
therefore be seen in a spectrum as a broad-line (BL) feature offset 
from a narrow-line (NL) feature, if the recoiling AGN is able to il- 
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Figure 13. Active lifetimes of recoiling BHs as a function of v^. Results are shown for the same six galaxy models shown in Fig. 1 1 1 [ separated by mass 
ratio for clarity (q = 1 on the left, q = 0.5 on the right). On the left: qlfg0.6a (red), qlfg0.3a (blue), qlfgO.la (green); on the right: q0.5fg0.6a (magenta), 
q().5fg0.3a (cyan), and q0.5fg0. la (gold). For each model, simulations were run for v-^lv csc = 0.4 — 1.2. In all plots, the solid and dashed lines denote the 
three different definitions of *agn discussed in the text, and the shaded regions denote the range of values between these bounds. The upper and lower solid 
lines denote *agn for Lbol > 3% and 10% -^Edd. respectively. The dashed lines show *agn for L^oi > L m in, where L m i n = 3.3 X 10 9 L . Unlike 
Fig. 1 1 1 1 here we have applied cuts to 4agn sucn that only the active time <agn < ^settle i s shown. In the top two windows, Lbol is calculated from the 
Bondi accretion rate only, while the bottom two windows show *agn with the ejected-disk accretion model included as well. 



luminate NL clouds as it leaves the center of the galaxy. Otherwise, 
such a system might have a spectrum with no NLs and BLs that 
are offset from the redshift of the host galaxy's stellar light. Be- 
cause the physics of the NL region are not included in our analysis, 
we consider the velocity offset between the recoiling BH and the 
the stellar center of mass. We calculate the AGN lifetime as before, 
based on the BH luminosity at each timestep, but in this case we are 
only interested in the lifetime of offset AGN; that is, we calculate 
the time for which the BH is active and exceeds a minimum veloc- 
ity offset AWnjn. We note that the limit Au m i n is imposed on the 
physical velocity difference between the BH and the host galaxy. 
The line-of-sight velocity seen by a typical observer will be lower, 
so the velocity-limited AGN lifetimes shown here are upper limits. 

Fig. [14] shows the results of our analysis for Av m i n = 
500 km s~ (top windows) and 800 km s~ x (middle windows). 
For Au m in = 500 km s _1 , the offset AGN lifetime for low- 
/ g as, low-Uk mergers is much smaller than the total AGN life- 
time in Fig. [T3] v csc is relatively low in these models, so kick 
speeds are also lower for fixed Vk/v asc . In the q0.5fg0. la merger, 
for example, v CBC = 777 km s _1 , so the BH never exceeds the 
500 km s _1 limit for v^/v CBC < 0.6. The offset AGN lifetimes for 



low fk/t'csc in gas-rich models are a factor of 5- 1 lower than the to- 
tal AGN lifetimes, but the peak lifetimes are still ~ 10 — 50 Myr. At 
high Uk/fcsc (> 0.8), there is little difference between t agn with 
or without the velocity cut. In these cases, the ejected disk domi- 
nates the BH accretion, and the BH speed is above 500 km s _1 for 
nearly all of the AGN phase. 

When Avmi n is increased to 800 km s _1 , the AGN lifetimes 
drop sharply for most merger models and kick speeds (see the mid- 
dle windows in Fig.|14), Only in the equal-mass, gas-rich models 
shown here (qlfg0.6a & qlfg0.3a) do the BHs appear as bright off- 
set AGN, and for only a narrow range of kick speeds. These two 
models have v csc > 1100, such that BHs kicked with ~ 800 km 
s _1 may have tightly-bound orbits. In these cases the BHs make 
many rapid pericentric passages through the dense central region, 
where they encounter an ample gas supply from which to accrete. 
Interestingly, in the qlfg0.6a and qlfg0.3a models, kicks < 800 km 
s _1 (< 0.6 v CS c) have nonzero lifetimes in this plot. These are ex- 
amples of the phenomenon discussed in § |4.2| in which the BHs are 
kicked before v eac reaches its maximum, so the BHs receive mild 
velocity boosts on subsequent passages through the increasingly 
deep central potential well. 
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Figure 14. The BH active lifetime is plotted in the same manner as in Fig. [13] and for the same galaxy models, but after minimum- velocity or minimum- 
separation cuts have been imposed. In all cases, both the Bondi accretion and the ejected-disk model are included in calculating Lbol- The top two windows 
show the timescale for which the BHs are active and also have a velocity greater than 500 km s —1 relative to the host galaxy's stellar center of mass. The 
middle two windows show *agn f° r a minimum velocity of 800 km s~ 1 . The bottom two windows show *agn f° r when the BH is spatially offset from the 
stellar center of mass by > 1 kpc. 



In the other (lower-g, / gas ) models, the central region is less 
dense, such that bound BHs have lower kick speeds and lower 
accretion rates during pericentric passages. In these models, only 
ejected-disk accretion is efficient at producing an offset AGN with 
Av > 800 km s _1 , which occurs only when the BHs are ejected 
entirely from the host galaxy. These ejected AGN have lifetimes < 
10 Myr and low luminosities (4agn < 1 Myr for L > 10% Z/Edd)- 



To be observed as a spatially-offset recoiling AGN, a recoiling 
BH must travel far enough from the galactic center to be resolved 
as distinct point source before exhausting its fuel supply. As with 
the velocity offsets, we define the spatial offsets relative to the stel- 
lar center of mass. The maximum galactocentric distance achieved 
by a recoiling AGN in any of our fiducial-mass simulations is ~ 15 
kpc by our lowest-luminosity AGN definition and only ~ 1 kpc by 
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the strictest definition. In our high-mass simulations, which have 
larger ejected disks, offsets > 50 — 100 kpc may be achieved while 
L > L m j n , but the maximum offset while L > 3 — 10% LEdd 
is still only a few kpc. Offsets of ~ 1 kpc correspond to an an- 
gular separation of ~ 0.55" at z = 0.1; these could be resolved 
with HST/JWST as long as they are not observed edge-on. Offsets 

> 3 — 4 kpc could be resolved at this redshift with SDSS and Chan- 
dra. To determine the lifetimes of the spatially-offset AGN in our 
simulations, we impose a minimum-distance cut of A_R m i n = 1 
kpc from the stellar center-of-mass, and we identify an offset AGN 
as one that simultaneously has AR > Ai? m i n and meets our AGN- 
luminosity criteria. The bottom windows of Fig.[l4]show the results 
of this analysis. We find that above Vk/fesc = 0.5 — 0.7 (physi- 
cally, the kick speed above which i? ma x exceeds 1 kpc), the offset 
AGN lifetimes can be long, but only for low-luminosity AGN. In 
these models, 4agn ranges from 2—110 Myr for L > L mm , and 
0.3 — 8 Myr for L > 3% I/Edd- In none of our simulations does 
the offset exceed 7? m i n for more than 1 Myr while its luminosity is 

> 10% I/Edd • The reason that spatially-offset AGN are only seen at 
low luminosities is as follows. For recoiling BHs with short-period 
oscillations, Bondi accretion is dominant, and most of the bright 
AGN activity occurs during pericentric passages with small AR. 
For BHs on long-period orbits or escaping BHs, the ejected disk 
accretion dominates. In this case the brightest AGN phase occurs as 
the BH leaves of the central dense region, and once the BH passes 
A_R m i n , its luminosity has already begun a monotonic decay. 

We can conclude from Fig.[l4]that in most of our models, low- 
luminosity recoiling AGN may be distinguishable via spatial off- 
sets > 1 kpc for kick speeds > 0.5 — 0.7 v csc , but as kinematically- 
offset AGN with A« m j n = 800 km s _1 , they may be distinguish- 
able only for kick speeds > w csc . These offset AGN generally 
have luminosities < 3%1/Edd, and their luminosity owes mainly 
to accretion from the disk ejected with the BH. The exceptions 
are equal-mass, gas-rich mergers (here, qlfg0.6a & qlfg0.3a), in 
which recoiling BHs on bound trajectories may experience multi- 
ple phases as bright (L > 10%LEdd), kinematically-offset AGN 
accreting from ambient gas during pericentric passages. This sce- 
nario can only occur in merger remnants with large central escape 
speeds ( > 1100 — 1300 km s" 1 ), such that BHs may receive kicks 

> 800 km s~ x and still undergo short-period oscillations. In lower- 
q, lower- /gas mergers, the central gas density is lower, and the kick 
speeds for bound trajectories are smaller, so this bright, velocity- 
offset phase does not occur. If kinematic offsets can be observed 
with a resolution of 500 km s _1 instead of 800 km s , velocity - 
offset-AGN lifetimes increase significantly, and measurable veloc- 
ity offsets can be attained for a much wider range of kick speeds. 

In the lowest q, / gas merger for which we have varied (the 
q0.5fg0. la model), the offset AGN lifetimes are always < 10 Myr 
for L > Lnm and < 2 Myr for L > 3%L E dd- Our q0.5fg0.1a 
model has 10% gas initially, but by the time of the recoil kick only 
~ 4% of the baryonic mass is in gas. We can conclude that merger 
remnants with gas fractions lower than this and with q<0.5 are 
unlikely to fuel bright offset AGN according to our criteria, at least 
for recoiling BHs with A/bh ~ 10 7 Mq. Larger BHs will have 
larger gas reservoirs and longer AGN lifetimes, and the reverse will 
be true of smaller BHs. For example, in the q0.5fg0. la model, a BH 
with Uk/^csc = 0.9 has a velocity-offset (Au m i n = 800 km s _1 ) 
AGN lifetime < 1 Myr by all AGN criteria, while in the 10 x more 
massive q0.5fg0.1M10x model, the velocity-offset ^agn is 27 Myr 
(for L > L m in). In general, larger BHs will have longer offset- 
AGN lifetimes for fixed Uk/v esc - However, larger BHs are found 



in larger galaxies, where a higher kick speed is required to reach 
an appreciable fraction of the escape speed, and such kicks have a 
lower probability of occurring. We discuss recoil kick probabilities 
further in § |6.1| 

As a final caveat, we note that when Bondi accretion is domi- 
nant, the offset- AGN lifetimes are subject to the assumption that no 
significant delay occurs between when material becomes bound to 
the BH and when it is actually accreted. The brightest AGN phases 
of this type occur during pericentric passages, when the velocity is 
high and the BH is subsequently traveling away from the galactic 
center. Therefore, substantial delays in gas accretion would likely 
decrease lifetimes of velocity-offset AGN and increase lifetimes of 
spatially-offset AGN. 

5.4 Mbh - <?* Relation 

Several authors have suggested that GW recoil may contribute to 
scatter in the observed black hole - host galaxy bulge relations 
due to ejected BHs ( Volonteri 20 07} or bound recoil ing BHs that 
leave the central galactic region l |Blecha & L oeb 2008 , Sijac kTet al.| 
[2^1^ . [\OMrtenl i |2007) used semi-analytic models to estimate the 
effect of ejected BHs at high rcdshifts on the BH mass - bulge stel- 
lar velocity dispersion relation at z = 0. In this work, we have al- 
ready demonstrated that, even when they remain bound to the host 
galaxy, recoiling BHs can undergo substantially different accretion 
histories than their stationary counterparts, and that the resulting 
AGN lightcurves may also vary significantly as a function of kick 
speed (see Figs. [10] [T3] [T4| >. Here, we attempt to quantify the ef- 
fect of individual GW recoil events on the Mbh — o"* relation. Be- 
cause we have simulated all of our merger models with both with 
no recoil kick and with = 0.9 v CBC , we are able to compare the 
Mbh — f* relations that result from each of these samples. The fi- 
nal BH masses and LOS-averaged stellar velocity dispersions from 
both sets of simulations are plotted in Fig. |15| the right window 
also shows simulations with v^/v csc 0.4 — 0.8. We include only 
models with the fiducial merger orientation angles (orbits a & w - 
z) in this plot to avoid heavily weighting the fit toward the (q, / gas ) 
combinations for which we have varied the orbit. 

By examining the mass deficits of these recoiling BHs rela- 
tive to the stationary BHs in our simulations, we can quantify the 
interruption to BH growth caused by GW recoil. Note that we are 
probing only the effect of GW recoil on the growth of kicked BHs; 
the actual BH mass deficit resulting from BHs merging in a cos- 
mological framework may differ from these values. In particular, 
the Vk = 0.9 v esc simulations result in BH wandering times of at 
least several Gyr, so in reality these BHs will be lost to the galaxy 
and never contribute to the Mbh — o~* relation. If the galaxy expe- 
riences subsequent mergers, especially minor mergers, the kicked 
BH may be replaced with a smaller BH from the incoming galaxy, 
creating a larger mass deficit. We also stress that the data plotted 
here are not meant to represent a random sample, as we have varied 
parameters in a systematic way. Accordingly, the relations shown 
are not expected to reproduce the observed Mbh — c* relation. In 
particular, the scatter in the ttk — relation (0.13 dex) is smaller 
than the observed scatter (~ 0.25 — 0.3 dex, e.g., Tremaine et al.| 



|2lMl|Aller & Richstone|2007| [Hopkins et al.|2007b| >, owing toa 
number of factors such as the single merger orbital configuration 
used in this plot, the constant feedback efficiency assumed for all 
simulations, and small range of total galaxy mass spanned by the 
majority of our models. Furthermore, although scaling each kick to 
0.9 Uesc essentially allows us to capture the maximum effect of a 
single recoil event on BH growth, this does not reflect the fact that 
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Figure 15. Left panel: A/bh — CT * relation for BHs with no recoil kicks (black circles) and with v-^lv aBC = 0.9 (red triangles). Mbh is in a U cases the BH 
mass at the end of the simulation (t mrg + 2.9 Gyr) and <r» is the stellar velocity dispersion averaged over 100 random sight lines, where error bars give the 
range of sampled values. All merger models with orbits a & w - z are shown. The solid black line is a least-square fit to the no-recoil data and the red dashed 
line is fit to the high-recoil data; the fit parameters are indicated on the plot. Right panel: same data as left panel, but also including data for the six merger 
models in which we have varied v^/v CBC . The blue diamonds show these results for v^/v CBC = 0.4 — 0.8, and the blue dot-dashed line is a fit to these data. 
The green triple-dot-dashed line is a fit to all data shown. 



in galaxies with lower escape speeds, large v^/v CBC recoil events 
will occur more frequently than in massive galaxies. Such a trend 
with galaxy mass would steepen the slope of the Mbh — cr* re- 
lation, though this would also depend on the distribution of other 
merger parameters, such as q and / gas , with varying galaxy mass. 
Despite these limitations in comparing our simulations to the ob- 
served Mbh — c* relation, our results provide substantial insight 
into how GW recoil affects the growth of kicked BHs and the inher- 
ent complexities of this process. The relative differences in Fig.|15| 
between the simulations with and without recoil are robust and can 
easily be understood in physical terms. 

It can be seen by visual inspection of Fig.[l5]that the effect of 
high-velocity recoils is an overall downward offset in normaliza- 
tion and an increase in scatter. Because the LOS-averaged values 
of a, do not change significantly in the presence of a recoil kick, 
we focus our discussion on variation in the final BH mass. 

The downward shift in Fig.[T5]reflects the fact that GW recoil 
always reduces the final mass of a merged BH. We quantify the 
"mass deficit" of recoiling BHs relative to stationary BHs by cal- 
culating the fractional mass difference for each recoil simulation 
relative to its = counterpart, (A/fl n (0) - Mfi n (Dk)) / Mfl n (0), 
where Mfl n is the final BH mass in each case. Fig. [16] shows the 
mass deficit plotted versus «k for the six merger models discussed 
previously. The deficits generally increase with v^, though the de- 
tails vary substantially between models. In the q0.5fg0.1a model, 
recoiling BHs with «k = v CBC have ~ 1/3 lower final mass than 
do stationary BHs, while in the qlfg0.3a model, a — v CBC recoil 
results in a BH almost five times smaller than its stationary counter- 
part. Note also that the curves in Fig. [T6] all level out at v^/v CBC — 
0.7 — 0.9; this marks the kick speed in each model above which 
recoiling BHs are unable to accrete more gas after their ejected 
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Figure 16. Fractional mass deficit of recoiling BHs versus stationary BHs 
as a function of v^/v eBc . Mfi n (0) and Affl n (tJk) are the final BH masses 
(at t mrg + 2.9 Gyr) for stationary and recoiling BHs, respectively. The plot 
legend lists the models shown. 



disk has been exhausted. The models with lower / gas achieve 
this limiting value at lower fk/^csc, owing to the lower gas drag 
and shallower potentials in these remnants. The blue diamonds in 
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Fig. [13] mark the positions of these varying-Uk simulations on the 
Mbh — cr* relation. The points corresponding to lower Vk/v esc lie 
close to the no-recoil (black) line; those with higher «k have lower 
BH masses. 

In Fig. [T5] the normalization offset between the = BH 
population and the Vk/Uesc — 0.9 population is 0.44 dex, and the 
offset between — and the total population is 0.22 dex. How- 
ever, it has been demonstrated that the normalization of the black 
hole - galaxy bulge correlations depends on the efficiency of AGN 
feedb ack <|Di Matteo et al.|2005| |Robertson et al.|2006c[ |Hopkins| 
[eTaLl[2007ar In our simulations, we assume a constant fraction 
(5%) of the BH's bolometric luminosity couples to the surrounding 
gas as thermal energy. Stationary accreting BHs with lower feed- 
back efficiencies would grow larger before heating the surrounding 
gas enough to slow or halt accretion. Depending on how much ac- 
cretion occurs prior to versus after the BH merger, higher feedback 
efficiencies could reduce the normalization offset between our re- 
coiling and stationary BH populations. In general, however, station- 
ary BHs will accrete more after the BH merger than recoiling BHs, 
so there is likely to be a non-zero downward shift in Mbh — 
normalization caused by GW recoil. 

As previously stated, the disproportionate effect of recoil 
events on galaxies with small v csc could have some effect on 
the Mbh — o"* slope. In our simulations, where we have scaled 
kick speeds to v csc , there is no significant difference between the 
Mbh — c* slope for recoiling versus stationary BHs. This is to 
be expected, because the depth of the central potential well scales 
more or less self-similarly with total galaxy mass. We do, however, 
see a trend toward larger BH mass deficits for higher-q mergers, 
which arises because gas is driven more efficiently to the galaxy 
centers during coalescence, and recoiling BHs "miss" a larger ac- 
cretion phase when they are kicked out of this central dense region. 

The variation in BH mass deficits between merger models re- 
sults in larger intrinsic Mbh — cr* scatter for BH populations that 
have undergone recoils. The correlation for the no-recoil sample 
has a scatter of 0.13 dex, and the v^/v csc = 0.9 correlation has a 
scatter of 0.24 dex. In other words, the high-recoil sample has al- 
most twice as much intrinsic scatter as the no-recoil sample. This 
increase in scatter is a universal consequence of GW recoil; recoil 
events essentially add an extra variable to the determination of fi- 
nal BH masses. As opposed to stationary BHs, the final mass of 
recoiling BHs depends not just on how much gas is available to be 
accreted, but also when this accretion occurs relative to the time 
of the kick. In mergers with no recoil, Mbh is a strong function 
of the amount of gas driven into the central region (i.e., the depth 
of the central potential well). The details and timing of the cold 
gas flow in galaxy mergers depend nontrivially on factors such as 
q, /gas, the galaxy orbits, and the star formation rate, so it is natural 
that larger intrinsic scatter is found in the final masses of recoiling 
BHs. If the kick speed is lower, such that the BH can settle back to 
the center while an ample gas supply is still present, the BH mass 
deficit will be reduced, but the details of this late-phase accretion 
add yet another element of unpredictability to the final value of the 
BH mass. 

In addition to the larger intrinsic Mbh — cr, scatter for a BH 
population with fixed v^/v csc , we also anticipate an increase in 
overall scatter caused by the range of kick speeds that would occur 
in a realistic population of BHs and by BHs that have undergone 
multiple recoils through their formation history. The contribution 
of GW recoil to the total Mbh — cr* scatter is difficult to predict 
from our simulations; quantitative predictions would require fol- 
lowing detailed BH formation histories in a cosmological frame- 



work, which is beyond the scope of this paper. However, we at least 
gain a sense of how a range of kick speeds affects the Mbh — cr, re- 
lation by plotting our varied-Uk simulations along with the i>k = 
and Wk/fcsc = 0.9 populations. The right panel in Fig.[T5jshows the 
same data as the left panel but also includes data for the simulations 
in which we have varied v^/v csc (blue diamonds). The blue line is 
a fit to these points, and the green line is a fit to all data shown. 
The fit to all the data has a scatter of 0.25 dex, which as expected is 
larger than the scatter for the Uk = population. 

The most robust of our conclusions from this analysis is the in- 
crease in intrinsic Mbh — &* scatter for BHs that have undergone 
a single recoil event, by a factor of ~ 2 in our simulations. This 
finding suggests that GW recoil events may be a non-negligible 
contribution to the scatter in the observed Mbh — cr„ relation, es- 
pecially because BHs at z = may have undergone multiple recoil 
events throughout their formation history. 

5.5 Evolution of Central Galactic Region 

Fig. [T7] compares the SFRs for simulations with no GW recoil 
(black dashed curve), with Vk/v csc = 0.9 (red solid curve), and with 
no BHs (blue dotted curve) in a gas-rich, equal-mass merger model 
(qlfg0.4a). It is clear that after the merger, the simulation with a 
large recoil kick has a higher SFR than the simulation with — 0. 
The post-merger SFRs are similar in the high-recoil case and in the 
simulation with no BH at all. This same behavior can be seen in 
the high-recoil simulations of model q0.5fg0.3a shown in Fig. |10| 
though to a lesser degree. These enhanced SFRs occur because 
when the BH is removed from the center of the merger remnant, 
the AGN feedback is also displaced. The lack of central energy 
input allows star formation to continue unhindered in this region 
for a longer time. We can actually see the effect of the recoiling 
AGN feedback on the instantaneous SFR in Fig.[T7] The SFR curve 
for the Vk/v csc — 0.9 simulation has a vaguely periodic shape; in 
fact, each small trough in the SFR coincides with a pericentric pas- 
sage of the BH through the galactic center, which imparts a small 
amount of feedback energy to the dense region. This supports the 
connection between the recoiling BH dynamics and the galactic 
SFR, but the main consequence of the recoil event is that star for- 
mation is never as strongly quenched as when the BH does not 
recoil. In the example shown in Fig. [T7] the merger remnant with 
a stationary BH has almost twice the amount of gas at the end of 
the simulation as the merger with v^/v csc = 0.9. This corresponds 
to an increase in total star formation of 3.4 x 1O 9 M0 in the recoil 
simulation, or ~ 3% increase in total stellar mass, solely due to the 
displacement of the BH and its feedback. 

In Fig. [T8] we can clearly see the effect of this enhanced star 
formation on the remnant stellar density profile. The figure shows 
the remnant stellar density profiles for the same three simulations 
shown in Fig. [IT] with «k = (black curves), Vk/v esc ~ 0.9 (red 
curves), and without BHs (blue curves). The results are further bro- 
ken down into the contribution from stars formed prior to the galaxy 
merger ("old" stars, thin curves) and those formed during merger- 
driven starbursts ("new" stars, thick curves). In all cases, the "new" 
stellar population dominates in the central kpc, but in the remnant 
with a recoiling BH the central stellar density is almost twice that 
of the remnant with a stationary BH. In fact, this stellar profile very 
closely resembles that of a remnant without a BH. Therefore, re- 
coil events in gas-rich mergers can prolong starburst phases, creat- 
ing denser and more massive stellar cusps. We have also examined 
the (unattenuated) U — B color evolution of the stellar component 
in these no-recoil, large-recoil, and no-BH simulations. The colors 
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in all three cases are similar; by the end of the simulation all have 
U — B ~ 0, i.e., these merger remnants are entering the "green val- 
ley" in transition from blue to red. However, in the v^/v csc = 0.9 
and no-BH simulations, which have prolonged star formation, the 
stars are slightly bluer than in the no-recoil simulation. While this 
difference is small, ~ 0.03 mag, it further demonstrates that in 
merger remnants with rapidly-recoiling BHs, the (blue) starburst 
phase may be prolonged, and hence the transition to red elliptical 
slightly delayed, solely due to the effect of GW recoil. In § |6.5[ we 
discuss the implications of these findings for observable properties 
of ultraluminous infrared galaxies (ULIRGs) and elliptical galax- 
ies. 

We note that this is the opposite effect of one that has been 
proposed for gas-poor mergers, where repeated BH oscillations 
through the central region of the merger remnant may scour out 
a stellar core {Boylan-Kolc hin et al.1 12004, Gualandris & Merritt| 
2008 1. Such core-scouring is expected to occur on scales of < 10s 



ot pc, which are below our resolution limit. However, we do com- 
ment that when even a small amount of gas is present in our merger 
simulations, it is efficiently driven to the center of the merger rem- 
nant. 



Figure 17. Total star formation rate (SFR) vs. time for the qlfg0.4a merger 
model. The black dashed line represents the = simulation, the red 
solid line represents the v^/v CBC = 0.9 simulation, and the blue dotted line 
represents a simulation with no BHs. The vertical line marks the time of BH 
merger in the simulations with BHs. 
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Figure 18. Final stellar density profiles (2.9 Gyr after t mr g) for the same 
three simulations shown in Fig. |17| with the same color coding. The thick 
curves are the "new" stellar populations that formed during the merger and 
that dominate within the central ~ kpc, and the thin curves are the "old" 
stellar populations from the initial stellar disks of the progenitor galaxies. 
In each case, the black curves are the profiles for the = simulation, 
the red curves represent the v^/v CBC = 0.9 simulation, and the blue curves 
represent the simulation with no BHs. Poisson error bars are shown on the 
densities in each radius bin. The vertical dotted line denotes the gravita- 
tional softening length (80 pc) in these simulations. 



6 DISCUSSION 

Using the SPH/N-body code GADGET-3, we have generated a suite 
of galaxy merger simulations both with and without a recoil kick 
applied to the central SMBH at the time of BH merger. Owing to 
the large range of parameters sampled, we were able to analyze 
systematic trends in recoiling BH behavior with variation in galaxy 
mass ratio, total galaxy mass, initial gas fraction, orbital configu- 
ration, and BH merger time. In addition, we have used a range of 
recoil kick velocities, which are scaled to the galactic central es- 
cape speed in each case. We followed the trajectories and accretion 
history of the recoiling BHs, as well as evolution of host galaxy 
properties such as the star formation rate and depth of the central 
potential well. For the BH accretion, the Bondi-Hoyle accretion 
rates were used, but we also included a time-dependent model for 
accretion from a gas disk carried along with the recoiling BH. 



6.1 Recoil Trajectories & Kick Probabilities 

The recoiling BH trajectories in our simulations are character- 
ized by low-angular-momentum orbits (i.e., "oscillations" about 
the galactic center) that are also highly three-dimensional. In other 
words, the baryonic component of the merger remnant dominates 
the BH orbits even if they extend well into the halo, but the unre- 
laxed nature of this remnant creates complicated BH trajectories. 
We find that the amplitude and duration of these oscillations varies 
widely between different galaxy models and kick velocities; some 
BHs settle back to equilibrium at the galactic center within a few 
Myr, while others may remain on large orbits through the halo for 
a Hubble time. 

The amplitude and duration of BH oscillations in a given re- 
coil event clearly depend on the central escape speed, i.e., the depth 
of the central potential well. However, by normalizing kick speeds 
in our simulations to v eac at the time of the kick, we are able to com- 
pare different galaxy models regardless of their escape speed. We 
find substantial variation in recoil trajectories between models even 
for fixed v^/v csc . Recoil oscillations are more readily suppressed in 
mergers with higher q and / gas . Such merger remnants are more 
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centrally concentrated and gas-rich, resulting in steeper central po- 
tential wells and increased gas drag and dynamical friction. 

In the models for which we have varied ^k/Wesc, we find that, 
with the exception of our collisionless simulations, recoiling BHs 
with «k/«csc = 0.4 — 0.7 reach galactocentric distances < 1 kpc 
and settle back to the center in < 0.5 Gyr. The corresponding range 
of kick speeds is 310-880 km s~ . 

For large recoil velocities, GW recoil events can produce BH 
oscillations that persist through the end of our simulations, 2.9 Gyr 
after the kick. The largest orbits can extend well into the galactic 
halo, with galactocentric distances up to ~ 100 kpc. In our models 
these long-lived oscillations occur for kicks > 0.6 — 0.9 v csc , again 
depending on the mass ratio and gas content of the merger. 

Furthermore, as illustrated in § |4.2| in gas-rich mergers with 
q ~ 1, the central escape speed may increase rapidly around the 
time of BH binary coalescence. In this case, the trajectory of the 
recoiling BH may also depend on the time of the BH merger rela- 
tive to the formation time of the central density cusp. Of all the fac- 
tors influencing the BH trajectories, this is the most unpredictable, 
as it depends on the relative timing of two short-timescale pro- 
cesses. This adds an element of uncertainty to any attempt, includ- 
ing the current study, to model populations of recoiling BHs in gas- 



rich galaxies. Our assumption of rapid BH mergers (t n 



^coal) 



means that our BH oscillation amplitudes may be upper limits in 
mergers with strong v csc evolution. However, the effect of this un- 
certainty on the results presented here is limited, due to the rela- 
tively small number of merger models affected and our use of the 
normalization v^/veac- 

We have also tested whether recoil trajectories are affected by 
the direction of the kick. Our generic merger remnants (those that 
are not formed from highly-aligned galaxy orbits) are highly dis- 
turbed at the time of the BH recoil and therefore do not have disk- 
like structures, i.e., there is no "special" direction in which a recoil- 
ing BH would clearly experience much greater drag. Even with the 
clumpy, irregular density profiles of our merger remnants, we find 
that the trajectories of recoiling BHs have little dependence on the 
kick direction. Much more important in determining the amplitude 
and duration of recoiling BH oscillations is the magnitude of the 
kick velocity relative to the escape speed, along with the aforemen- 
tioned global galaxy properties such as the initial gas fraction and 
mass ratio. We note that there should indeed be a dependence of 
the recoil trajectory on the kick direction with respect to the plane 
of the accretion disk; recoils directed into the disk plane likely ex- 
perience more gas drag and have shorter oscillation timescales. In 
addition, (Rossi et al. 2010) have shown that the ejected disk mass 
and AGN luminosity may depend strongly on the direction of the 
recoil kick with respect to the accretion disk plane. We do not re- 
solve this small-scale region in our simulations; however, recoil 
kicks > 500 km s _1 will always be nearly perpendicular to the or- 
bital plane of the binary ( Campa nelli et al.|2 007a b ; Lo usto & Zlo-| 
|chower|2009[ l, which is unlikely to be perpendicular to the orbital 
plane of the accretion disk. Therefore, for the recoil velocities we 
consider ( > 300 km s -1 ), we can assume that the BH is unlikely 
to be kicked directly into the accretion disk plane. 

We mention in § [2] that according to the results of NR simu- 
lations, the recoil kick speeds in our fiducial-mass merger models, 
which we have scaled to the escape speeds of ~ 700 — 1250 km 
s" 1 , should occur with fairly high probabilities. Even if the BH 
spin magnitudes are randomly distributed between and 1, 41.6% 
of major mergers should have kicks above 50 km s -1 , and 13.5% 

(2010b) 
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should be above 1000 km s . Furthermore, 
have conducted statistical analysis of GW recoil spins and kick 



speeds using post-Newtonian approximations calibrated to numeri- 
cal relativity results. They provide kick probability distributions for 
several mass ratio bins, which we can use to estimate the relative 
probability of "high" (v^ ~ v csc ) or "moderate" (v^ ~ 0.6 iw) 
kick speeds in various merger models. Coincidentally, for mass ra- 
tios 0.3 < q < 1.0 and 300 < v k < 1000 km s"\ which rep- 
resent a large fraction of our fiducial-mass simulations, the kick 
probability distribution is quite flat. For these models, the relative 
probability of a given kick speed between models, as well as the 
relative probability of high versus moderate v k in a given model, is 
roughly unity. In our highest-mass (qlfg0.3M20x) simulation, the 
probability of a kick near v CEC — 2500 km s _1 is about 10 times 
smaller than the probability of a kick with 0.6 v csc and is also ~ 10 
times smaller than the kick probability for Uk ~ v csc in the qlfg0.3a 
model. For mergers with q < 0.25, which we do not model in our 
study, the kick distribution becomes heavily skewed toward veloc- 
ities < 500 km s" 1 . 

An important caveat to the interpretation of recoil kick prob- 
abilities is that the relevant velocity for observations is the line- 
of-sight (LOS) velocity, not the actual velocity. For random lines 
of sight, the probability of observing a given kick i>k as a kine- 
matic offset > An m i n is therefore reduced by a factor P(Awlos > 
Aw m i n ) = 1 — A^min/uit. For example, about 20% of recoils 
with D k = 1000 km s" 1 will have Aulos > 800 km s _1 . 
|Lousto et al.|(2010b| l find that the overall probability distribution 
for Adlos is lower than the i>k distribution by a factor of between 
~ 2 (for kicks > 500 km s _1 ) and ~ 20 (for kicks > 2500 km 
s _1 ). Similarly, the probability of observing a recoiling BH with 
spatial offset Ai? LO s > Ai? min is yjl- AR 2 mh J AR 2 ; i.e., a 
majority of spatial offsets > 1.2i? m i n will have LOS offsets ex- 
ceeding the minimum observable offset. 



6.2 Recoiling AGN Lightcurves 

The final coalescence of two merging, gas-rich galaxies typi- 
cally triggers both a burst of star formation and a phase of rapid 
BH accretion. Assuming the BHs themselves coalesce on a short 
timescale, the resulting recoil event coincides with a period of high 
activity in the central region of the merger remnant. Thus, we ex- 
pect recoiling BHs in gas-rich mergers to have substantially dif- 
ferent accretion histories than their stationary counterparts. This is 
indeed the behavior seen in our simulations, with greater disruption 
to the BH accretion history for higher kick velocities. 

Recoil events that produce short-period BH oscillations 
(Vk/v e sc ~ 0.4—0.8, depending on the merger model) cause the BH 
accretion rate to drop at the time of the kick, thereby reducing the 
peak luminosity of the BH. However, the BH remains in the central 
few kpc of the merger remnant where gas densities are large enough 
that Bondi-Hoyle accretion is efficient. Because the BH's feedback 
energy is deposited in a much larger volume than for stationary 
BHs, it is unable to heat the gas sufficiently to completely termi- 
nate the AGN phase. The recoiling BH may therefore experience a 
longer active phase by a factor of a few than its stationary counter- 
part, although at lower luminosities. This phase continues in some 
cases until the end of the simulation, 2.9 Gyr after the BH merger. 
Owing to the lower luminosities of these extended AGN phases, we 
find that they do not increase the final BH mass relative to station- 
ary AGN. However, this type of accretion episode is distinct from 
other merger-triggered AGN phases in that it is not accompanied 
by a simultaneous starburst. The copious amounts of dust produced 
by starbursts can easily obscure an active central source, and the 
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starburst luminosity may overwhelm the AGN luminosity. Recoil- 
ing AGN that continue accreting long after the starburst is complete 
therefore may be more easily observable. 

When recoil kicks are large enough to produce long-period 
( > 100 Myr) BH oscillations (v k /v CBC > 0.6 — 0.9 depending on 
the merger model), the Bondi-Hoyle accretion rate drops precipi- 
tously as the BH is rapidly ejected to a lower-density region. Even 
though the BH remains on a low-angular-momentum orbit and re- 
turns several times to the dense central region, its velocity is gener- 
ally too high during these pericentric passages to allow for substan- 
tial gas accretion, and the gas is depleted by star formation between 
subsequent passages. 

In such high-velocity recoil events, a portion of the inner ac- 
cretion disk may remain bound to the BH. We account for this by 
assigning to each recoiling BH a time-dependent, isolated accretion 
disk with a mass and initial accretion rate based on the BH mass, 
recoil velocity and the accretion rate at the time of the kick. The 
mass of these ejected disks is only a few percent of the BH mass, 
so they do not contribute significantly to the final BH mass. How- 
ever, the accretion from these disks does greatly extend the active 
lifetime of rapidly-recoiling BHs in the regime where Bondi-Hoyle 
accretion is inefficient. We find that the AGN lifetime in our models 
may increase from <C 1 Myr to ~ 20 Myr when ejected-disk accre- 
tion is included, though these are generally low-luminosity AGN. 
The peak lifetimes in the Bondi-Hoyle regime are generally longer 
and correspond to brighter AGN, but for a fairly narrow range of 
kick speeds. £ agn in the Bondi-Hoyle regime may be up to ~ 300 
Myr for low-luminosity AGN and between ~ 10 — 100 Myr for 
AGN accreting at 3 - 10%L E dd- 

6.3 Offset AGN Lifetimes 

From an observational standpoint, our primary interest is not the to- 
tal active lifetime, but rather the time for which the BH is active and 
can be distinguished as a recoiling BH. We find that kinematically- 
offset AGN - those that are luminous enough to be detected and 
have velocity offsets large enough to be spectrally resolved - occur 
in two distinct physical scenarios. In the first scenario, the recoil- 
ing AGN travel far from the galactic center and are powered by 
accretion from the ejected disk. This occurs for recoil velocities 
near or above the central escape speed. In our simulations, AGN 
with velocity offsets > 800 km s _1 have lifetimes < 10 Myr in 
most cases. In most of our galaxy models, recoiling AGN can only 
meet this velocity criterion for > v csc . For reference, the SDSS 
quasar study conducted by Bonning et al. (2007), which yielded a 
null result, had a spectral offset limit of 800 km s _1 . 

The other physical mechanism for producing kinematically- 
offset recoiling AGN is via multiple pericentric passages through 
a central dense region. Note that marginally-bound recoiling BHs 
may also experience close passages, but these will be rapid and few 
in number; the cumulative probability of observing such an event 
is low. If a recoiling BH is able to settle back to the galaxy center 
within ~ 0.5 — 1 Gyr, it will undergo many short-period oscilla- 
tions through the central region. In this case, the BH may be able 
to accrete efficiently from the ambient gas during high-velocity 
close passages, producing a kinematically-offset AGN. This sce- 
nario favors massive, gas-rich galaxies with high central densities: 
the larger the central escape speed and the steeper the potential 
well, the higher the velocity with which a BH can be kicked and still 
remain tightly bound to the galaxy. In our q = 1, gas-rich mergers, 
offset AGN with Av > 800 km s _1 can have lifetimes up to 100 
Myr, and may have lifetimes > 10 Myr at bright (> 10% LBdd) 



luminosities. Especially for these AGN, an increase of spectral res- 
olution to achieve a limit of 500 km s _1 could make a substan- 
tial difference; in our models, the lifetimes increase by a factor of 
~ 5 — 10 with this lower value of A« m i„. 

The maximum spatial offset achieved by a recoiling AGN in 
our high-mass simulations is AR — 112 kpc by our most gen- 
erous AGN definition (L > 1.3 x 10 43 erg s" 1 ) and 3.2 kpc by 
our strictest definition (> 10%1/EddX but for our fiducial-mass 
simulations these maximum offsets are only 15 and 1.3 kpc, re- 
spectively. To calculate spatially-offset AGN lifetimes, we require 
a minimum offset Ai? m i n = 1 kpc. Such an offset could be re- 
solved by HST or JWST at z<0.1, or by SDSS or Chandra at 
z < 0.03. In our simulations, recoiling BHs achieve i? ma x > 1 
kpc for Vk/v csc > 0.5 — 0.7. AGN lifetimes may be quite long near 
this threshold kick speed (~ 100 Myr), though only at low lumi- 
nosities. At higher kick speeds, lifetimes are generally ~ 1 — 10 
Myr and slowly decrease with higher Wk, as less gas bound to the 
recoiling BHs. 

Based on our calculated AGN lifetimes, we conclude that for 
escaping BHs (v^ > v csc ), probabilities for observing recoils via 
spatial offsets and kinematic offsets are similar. For recoiling BHs 
on bound trajectories, velocity-offset AGN with high luminosities 
are favored for massive, gas-rich galaxies (q = 1, / gas = 0.3 — 0.6 
in our models). Spatially-offset AGN with lower luminosities have 
similar lifetimes for all models. Recoiling BHs in merger remnants 
that are both gas poor (/ gas < 0.1 initially) and result from lower- 
q mergers (q < 1/3) are not expected to produce long-lived offset 
AGN. 



6.4 M BH - Scatter & Offset 

We have demonstrated that recoiling BHs have lower final masses 
than their stationary counterparts, though the mass deficit varies 
considerably with kick speed and merger remnant properties. This 
may affect the formation of the observed correlations between 
SMBHs and the bulges of their host galaxies. In particular, we con- 
sider the possible effects of these mass deficits on the observed 
Mbh — c* relation by comparing the A/bh — er» correlations of our 
sets of «k = and «k = 0.9 v csc simulations, as well as our simula- 
tions with varying Wk. We do not find any evidence for a change in 
Mbh — o"* slope caused by recoils with fixed v^/v esc , but in a realis- 
tic population of galaxy mergers, galaxies with smaller v csc would 
be disproportionately affected by recoils with «k near or above v esc . 
This could contribute to some steepening of the Mbh — o* slope. 
In our simulations, we find that GW recoil events may contribute 
to a downward shift in normalization, an increase in overall scatter, 
and an increase in intrinsic scatter of the A/bh — a* relation (see 
Fig.|15}, The magnitude of these effects, particularly the normaliza- 
tion offset and total scatter, are sensitive to the galaxy population 
considered and the detailed merger and recoil histories of the cen- 
tral BHs. Therefore, we cannot make quantitative predictions about 
the contribution of GW recoil to the observed normalization and 
scatter of the A/bh — a* relation; this would in fact be difficult 
to predict with any method. A more realistic approach involving 
BH and galaxy populations derived from a cosmological frame- 
work would still be plagued with uncertainties in the recoiling BH 
masses. We have shown that these masses are sensitive numerous 
factors that do not lend themselves to semi-analytical modeling, 
including the (possibly evolving) central potential depth, the BH 
merger timescale, and the detailed BH accretion history. 

Despite these uncertainties, the existence of a recoiling-BH 
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contribution to the A/bh — cr* correlation and its scatter is a nec- 
essary consequence of GW recoil events. In particular, we demon- 
strate that a universal effect of GW recoil is to increase the intrinsic 
scatter in BH masses. The disruption to BH accretion caused by a 
recoil event introduces an extra variable into the determination of 
the final BH mass, namely, the timing of the recoil kick relative to 
the merger-induced burst of BH accretion (i.e. quasar phase). In our 
recoil simulations, the intrinsic Mbh — o"* scatter of our rapidly- 
recoiling BH population (vk/v csc — 0.9) is almost 2x the scatter 
in our stationary-BH population. This is the contribution to scatter 
caused solely by interruption of BH growth by a recoil kick; other 
factors not accounted for here could easily increase the contribution 
of GW recoil to A/bh — o"* scatter. While a more realistic popula- 
tion of BHs will have smaller average mass deficits resulting from a 
range of kick speeds, the BHs are also likely to have undergone nu- 
merous mergers and recoil events in the past, which would increase 
the overall scatter. Subsequent merger events, in which an incoming 
BH may replace an ejected BH, could also increase scatter. The rel- 
atively small scatter of the empirical relation, 0.25 - 0.3 dex, places 
constraints on the contribution of GW recoil and, hence, on the 
frequency of large v^/v CBC recoil events. Similarly, a non-zero con- 
tribution of GW recoil to the observed scatter places constraints on 
the contribution from other sources, such as redshift evolution of 
the M BH - a* relation (e.g.,|Haehnelt & Kauffmann|2000)|Shields] 
|et al.|2003||Robertson et al.|2006c| >. 

6.5 Broader Implications: Starbursts, ULIRGs, and 
Elliptical Galaxies 

Previously, a model for galaxy evolution has been outlined in which 
mergers play a central role (e.g., Hopkins et al. 2006a, 2008b 
|Somerville et al.||2008[ >. In this picture, galaxies form originally 
as disks, and the population of ellipticals is built up through time 
from mergers. Gas-rich mergers pass through intermediate phases 
in which they host significant star formation and black hole activity, 
first as ULIRGs, then as quasars, and finally as dormant elliptical 
galaxies. 

There is a great deal of observational support for this scenario. 
Nearby ULIRGs are invariably associated with gas-rich mergers 
(Sand erset al.| 19 88a; San ders & Mirabel|1996) and contain larg e 
quantities of molecular gas in their centers ^Sanders et al.|1991| >, 
driving nuclear starbursts. Many display evidence for central AGN 
through their "warm" spectra (Sanders et al. 198Kb). leading to the 
interpretation that quasars are the descendants of ULIRGs. As the 
starbursts and quasar activity fade, the merger remnant evolves into 
a passive elliptical (Toomre & Toomre 1972; Toomre 1977), mak- 
ing the transition from blue to red. The relics of this process can 
be seen in the central light concentrations in elliptical galaxies that 
were put into place during the starburst phase of their creation (e.g., 
Mihos & Hernquist 1994a, Hopkins et al. 2008c), and correlations 
between supermassive black holes and e.g. the central potentials of 
their hosts ( |Hopkins et al.|2007b||Aller & Richstone|2007| l. 

These fossil signatures, and the timing of the events that led to 
their formation, are sensitive to the interplay between star forma- 
tion and black hole growth and the feedback processes that accom- 
pany each of these. In principle, black hole recoil could impact all 
of these phenomena to varying degrees, providing constraints on 
the modeling we have presented here as well on the interpretation 
of the observations more generally. 

For example, some ULIRGs display evidence for significant 
AGN activity in the central starburst region. This is strong evi- 
dence that the BHs in these systems did not undergo rapid recoil 



events, and that any recoil motion was quickly damped out. How- 
ever, merger-triggered central AGN, particularly those in ULIRGs, 
may be dust-obscured or dwarfed by the luminous starburst. AGN 
that are instead kicked out of the central region by a large recoil 
may be more easily detected, provided a large enough region of 
the galaxy is observed. Additionally, we have shown that in some 
cases, bound recoiling AGN may have longer lifetimes than sta- 
tionary AGN. Such AGN may continue shining after the starburst 
is complete, which may also make them easier to observe. 

Moreover, in the simulations, at least, feedback from black 
hole growth plays a significant role in quenching star formation ( |Di| 
|Matteo et al.|2005 ; Springel et al. 2005a), allowing the remnant to 
evolve from a blue, actively star-forming galaxy, to a red, quiescent 
object. During this transition, the remnant can, for some period of 
time, be classified as a K+A or E+A galaxy l |Snyder|2010") . If the 
central BH were ejected, star formation would linger at higher rates 
than for a stationary BH, extending the duration of the blue to red 
transition, modifying the stellar populations, and ultimately yield- 
ing an elliptical with denser, more massive central stellar cusps. If 
the BH did not eventually return to the nucleus, subsequent gas- 
poor mergers would be less strongly affected by BH scouring, ob- 
scuring the relationship between core and cusp ellipticals of the 
same mass (Hopkins & Hernquist 2010b). 

Fig.[T8]shows that if the BHs recoil at sufficiently high speed, 
the remnant relaxes to a state similar to what would have happened 
had the galaxies not had black holes. Compared to the case of no 
BH recoil, the remnants with recoiling BHs (or no BHs) have more 
massive central stellar cusps. These more massive stellar cusps re- 
flect a more extended period of active star formation following fi- 
nal coalescence of the progenitor galaxies, lengthening the period 
of time to make the blue to red transition. Furthermore, the reduced 
efficiency of AGN feedback means that the central regions would 
remain more heavily dust-obscured than if there were no BH re- 
coil, making the remnant less likely to be visible as a K+A galaxy 
l |Snyder|2010) and potentially significantly reducing the number of 
K+A galaxies expected to be observable as merger remnants. 

Using our suite of merger models, we can make some qual- 
itative predictions about the morphologies of recoiling BH hosts. 
Fig. [3J shows that our merger remnants have some common mor- 
phological features. In general, at the time of BH merger (and re- 
coil) the merging galaxies have undergone final coalescence (i.e., 
they no longer have two distinct cores), but are still highly dis- 
turbed and have prominent tidal features. This is true even when 
the progenitor galaxies have low initial gas fractions (4%) or rela- 
tively small ratios (0.25). For minor mergers with q < 0.25, the dis- 
ruption to the primary galaxy would be less pronounced, but such 
mergers would also yield small GW recoil kicks. Highly-aligned 
(i.e. coplanar) galaxy merger orbits may result in remnants with 
regular disk structures (see Fig. 13), but these orbits should com- 
prise only a small fraction of major galaxy mergers. Therefore, we 
conclude that the majority of galaxy merger remnants that produce 
large GW recoil events have single cores but are still tidally dis- 
turbed at the time of the recoil kick. Furthermore, we have shown 
that the offset-AGN lifetimes for recoiling BHs in our simulations 
are < 100 Myr and are in most cases ~ 1 — 10 Myr. Tidal struc- 
tures in these major merger remnants generally persist for at least a 
few hundred Myr, so it follows that most observable recoiling AGN 
should be found in unrelaxed remnants that may still have visible 
rings or tidal tails. 
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6.6 Comparison to Recent Work 

Shortly before this paper was completed, Guedes et al. (2010) and 
|Sijacki et al.| ( |2010| > completed independent studies of recoiling 
BHs in gaseous mergers. However, both of these papers are ex- 
ploratory in nature, in that only a few examples of galaxy merg- 
ers are studied. While there is agreement between these studies on 
some fundamental conclusions, our large parameter study allows 
us to expand upon these findings, connect them into a more coher- 
ent picture, and present some entirely new results. Here we give a 
brief comparison of this work to the work of Guedes et al.| j2010) 
& |Sijacki et aT~| ( [2010l >. 

|Guedes et aTT| pOlO) use one equal-mass and two minor 
merger remnants as initial conditions for their simulations. The 
mass and spatial resolution used are comparable to ours, though 
in their equal-mass merger, higher resolution is obtained at late 
stages via particle splitting. The merger orbits were coplanar and 
prograde, a configuration that results in highly disk-like remnants 
(see our Fig.[5|. Recoils with a range of kick speeds were applied 
to the BH, and after a brief initial phase, the BH trajectories were 
followed semi-analytically. We note that they use v^/v csc as low as 
0.14, while we use only > 0.4 v eac , such that the BHs do not 
spend most of their time at radii ~ R BO tt- 

|Sijacki et al.| ( |2010"| > use isolated galaxy models as well as one 
galaxy merger model for their recoiling BH simulations. Like our 
simulations, theirs were performed using GADGET-3, with similar 
spatial resolution. However, they use only one merger model, which 
is massive and gas-rich, resulting in a remnant with a very high cen- 
tral escape speed (3500 km s _1 ) and steep central potential. This 
model is therefore fairly extreme, in that recoils with ~ v eac are 
required to eject the BH beyond the central few kpc. 

Despite the substantially different approaches of these two 
studies and our own, a common finding is that gas inflow during 
galaxy mergers can create centrally-dense remnants that impart sig- 
nificant gas drag to recoiling BHs, thereby reducing the amplitude 
and return time of their trajectories. That this effect remains impor- 
tant even for modest gas fractions ( < 0.1) illustrates the limitations 
of modeling GW recoil dynamics in purely collisionless systems. 

|Guedes etal . (2010) find that in their minor mergers, recoiling 
BHs exceed 600 km s _1 only for > v csc , and they do not cal- 
culate AGN lifetimes for escaping BHs. Kinematically-offset AGN 
have li fetimes up to a few Myr in their q — 1 merger. Sijacki et al. 
( |2010| note that velocity offsets up to 500 km s 1 may occur dur- 
ing the bright AGN phase in their isolated disk model, but do not 
discuss offset AGN in their merger model. 

Although they use similar values of L m i n and Afl m j n , |Guedes| 
|et al.| ([20 10) calculate longer spatially-offset AGN lifetimes than we 
find in our models. The disparity in offset AGN lifetimes may arise 
partly from the different dynamical treatments; Guedes et al.H2010) 
calculate the recoil trajectories semi-analytically, as opposed to our 
N-body approach. Further, they assume Mbh = M^dd for the 
ejected-disk accretion rather than a time-dependent, isolated-disk 
model. They argue that spatially-offset AGN are naturally longer- 
lived than kinematically-offset AGN because the former occur near 
apocenter and the latter near pericenter. However, we show that in 
some cases, kinematically-offset AGN may have longer lifetimes, 
because the BH is much more likely to be actively accreting during 
passages though the central dense region. 

As we have found in our recoil simulations, |Sijacki et al.| 
( |2010fr show that recoil events interrupt a phase of rapid BH ac- 
cretion in their gas-rich merger, causing a BH mass deficit relative 
to a stationary BH. They note that BH mass deficits may introduce 



scatter into BH mass - host galaxy relationships, which we demon- 
strate robustly in § |5.4| 

Furthermore, both our results and those of Sijacki et al. (2010) 
show that the central SFR is higher after a large recoil, due to the 
displacement of AGN feedback. Sijacki et al. (2010) note that this 
further depletes the gas supply available to the BH once it returns 
to the center. We demonstrate that SFR enhancement due to re- 
coil is strongest in gas-rich, nearly-equal-mass mergers, though it 
is present in all of our simulations to some degree, and we have 
outlined several possible observational consequences in § |6.5| 

In addition to connecting and expanding upon the results of 
|Guedes et allpOTo) and |Sijacki etal^pOTo) , many of our results 
are entirely new to this work. As we have discussed these previ- 
ously, we list them here only briefly. First, because we allow the 
BH merger time to be a free parameter, we are able to demon- 
strate that recoil trajectories may depend quite sensitively on the 
BH merger time in q ~ 1, gas-rich mergers. We eliminate much of 
this uncertainty in our simulations by scaling to u csc (tmrg). This 
approach also enables a much clearer comparison between models. 
Also unique to our work is that we test a range of recoil kick di- 
rections in generic merger remnants, demonstrating that in general, 
the kick direction is less important for recoil trajectories than the 
kick magnitude. 

Another novel feature of our models is the inclusion of a time- 
dependent accretion disk around the recoiling BH. We also calcu- 
late recoiling AGN lifetimes using three different luminosity lim- 
its, which allows us to differentiate between bright and faint AGN. 
We learn from this analysis that some recoiling BHs may actually 
have longer AGN lifetimes than stationary BHs, though at low lu- 
minosities. We also show that kinematically-offset AGN occur via 
two different mechanisms: repeated pericentric passages in mas- 
sive, gas-rich remnants, and ejected-disk accretion in high-velocity 
recoils (i>k > v csc in our models). We demonstrate that spatially- 
offset AGN can occur for a wide range of kick speeds and merger 
models. 

Finally, we use our suite of simulations to construct Mbh — 
relations for stationary and rapidly-recoiling BH populations. We 
find that GW recoil introduces a normalization offset and larger 
overall and intrinsic scatter into this correlation. 

We have shown that recoiling BHs may be observable as offset 
AGN; discoveries of these objects would both inform models of 
hierarchical galaxy formation and constrain event rates for future 
gravitational-wave observatories. In the meantime, the richly varied 
effects of GW recoil on SMBHs and their host galaxies are not to 
be discounted as important components of merger-driven galaxy 
evolution. 
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